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FOREWORD 

This  report  shows  tutorial  derivations  (a  form  of  student’s  guide)  of  all  key  equations  needed 
in  the  research  of  three-dimensional  eye  rotations,  which  happens  to  coincide  with  virtually  all  key 
equations  needed  to  understand  the  rotational  motions  of  spacecrafts  and  missiles.  Considering  the 
fact  that  mathematics  may  not  be  the  primary  field  of  competence  for  most  medical  researchers,  the 
author  of  this  report  took  special  care  to  explain  the  derivations  from  a  very  elementary  level  and  to 
gradually  progress  to  advanced  equations  in  a  logical,  self-containing  way.  The  author  introduces 
some  innovative  ways  of  explaining  otherwise  difficult  concepts  and  derivations  in  several  topics. 

A  few  years  back,  the  author  was  a  visiting  scientist  for  a  year  at  the  Department  of  Neurology 
of  the  Johns  Hopkins  School  of  Medicine  at  the  invitation  of  his  former  advisor,  Professor  David  A. 
Robinson.  During  his  stay  at  Hopkins,  he  realized,  based  on  requests  and  encouragement  from 
colleagues,  a  need  for  a  comprehensive  book  of  this  nature  containing  all  key  equations  and  their 
derivations  in  one  volume.  In  writing  this  report,  the  author  followed  the  pedagogical  philosophy  that 
the  best  way  to  understand  the  equations  and  to  be  able  to  use  them  with  confidence  is  to  be  able  to 
derive  them  from  the  basics  in  an  easy-to-follow,  yet  rigorous  way.  For  this  reason,  this  report  starts 
with  elementary  trigonometry  by  design,  and  advances  logically  to  a  higher  level. 

As  the  bibliography  provided  at  the  end  of  this  report  reveals,  all  equations  in  the  report 
originated  from  applications  to  spacecraft  and  missile  dynamics.  The  only  deviation  in  this  report 
is  that  the  Head  frame  is  identified  with  the  reference  frame,  while  the  Eye  frame  is  identified  with 
the  moving  or  rotating  frame. 

Although  this  report  starts  from  a  modest  level,  it  advances  toward  the  end  to  fairly  advanced 
esoteric  equations,  which  even  most  aerospace  engineers  may  not  be  aware  of,  or  may  not  encounter 
during  their  career.  These  equations  are  highly  desirable  for  a  deeper  insight  into  angular  rotations. 
For  this  reason,  although  this  report  is  written  primarily  with  medical  researchers  in  mind,  it  is  equally 
recommended  for  scientists  and  engineers  in  the  aerospace  field  dealing  with  the  rotational  dynamics 
of  missiles,  spacecraft,  and  aircraft  in  the  Department  of  the  Navy. 

Finally,  the  author  of  this  report  would  like  to  dedicate  this  modest  work  of  his  to  the  benefit  of 
medical  researchers  and  students  in  the  United  States  and  abroad  who  may  be  struggling  with 
unfamiliar  mathematics  in  their  work.  This  report  is  not  copyrighted  so  that  it  may  be  reproduced 
freely  for  educational  purposes. 

Questions  or  comments  concerning  the  contents  of  this  document  should  be  addressed  to 
NSWCDD,  Attn:  Dr.  Kee  Soon  Chun,  K44,  Dahlgren,  Virginia  22448-5100,  or  telefaxed  to 
540-653-8382,  or  e-mail  to  <kChun@nswc.navy.mil>. 


Approved  by: 

C.  A.  KALIVRETENOS,  Head 
Strategic  and  Strike  Systems  Department 
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PREFACE 


Dr.  Kee  Soon  Chun  was  first  exposed  to  the  oculomotor  system  around  1976  as  a  graduate 
student  in  my  laboratory  in  the  Wilmer  Eye  Institute  at  the  Johns  Hopkins  School  of  Medicine  where 
he  did  his  thesis  research  for  a  doctorate  in  Electrical  Engineering.  We  made  a  model  of  the 
vestibule-ocular  reflex  with  a  quick-phase  mechanism  to  generate  nystagmus  (A  model  of  Quick 
Phase  Generation  in  Vestibule-Ocular  Reflex,  by  K.  S.  Chun  and  D.  A.  Robinson,  Biological 
Cybernetics  28,  pp.  209-221(1978)).  It  introduced  the  concept  of  two  functions  of  time:  (1)  The 
“when”  curve  which,  when  the  eye  position  reached  it,  initiated  a  quick  phase;  and  (2)  the  “where” 
curve  to  which  the  quick  phase  carried  the  eye.  These  are  concepts  that  have  been  found  useful 
subsequently. 

Many  years  later,  in  1996,  Kee  Soon  was  granted  a  sabbatical  by  the  Navy.  After  working  for 
many  years  on  the  mathematics  underlying  the  flight  of  intercontinental  ballistic  missiles,  he  was 
motivated  by  a  desire  to  turn  his  talents  to  something  of  a  more  humanitarian  nature  and  thought 
again  of  Johns  Hopkins  and  our  laboratory.  He  joined  us  once  again  and  helped  our  work  on  a  neural 
network  model  of  the  neural  integrator  (the  one  that  converts  vestibular  velocity  commands  to 
oculomotor  position  signals).  While  doing  so,  he  became  aware  that  much  of  the  research  in  our 
laboratory  (that  of  Dr.  David  Zee)  involved  the  analysis  of  eye  movements  in  all  three-dimensions. 
This  subject  has  been  growing  in  importance  over  the  last  fifteen  years,  when  one-dimensional 
analyses  were  felt  to  be  largely  understood,  and  it  was  realized  how  easy  it  was  to  measure  torsion 
movements  with  the  eye-coil/magnetic-field  method.  This  has  lead  to  the  need  for  more 
sophisticated  mathematics  when  dealing  with  such  things  as  coordinate  transformations.  Listing’s 
Law,  and  quaternions.  All  of  this  caught  Kee  Soon’s  fancy  because  it  was  not  so  very  different  from 
what  one  needs  in  studying  the  flight  of  a  missile  over  a  revolving  planet.  He  was  inspired  to  bring 
these  mathematical  tools  together  and  catalogue  them  for  the  convenience  of  others.  Thus,  this 
publication  came  about. 

The  result  might  be  regarded  as  a  reference  source;  a  place  to  go  and  remind  oneself  of  the 
mathematical  bases  of  the  various  tools  we  use  in  oculomotor  analyses.  For  those  not  already 
comfortable  with  vectors  and  matrices,  this  text  may  be  too  terse  and  compact  to  serve  as  a  learning 
tool  although  the  sections  on  the  Fick’s  and  Helmholtz’s  coordinate  systems  may  help  the  beginning 
student  to  understand  these  potentially  confusing  representations.  Also  the  description  of 
semicircular  dynamics  (Appendix  B)  is  a  useful  introduction  to  anyone  beginning  a  study  of  the 
vestibule-ocular  reflex. 


David  A.  Robinson,  Professor 
School  of  Medicine 
Johns  Hopkins  University 

January  1999 
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SECTION  1 
INTRODUCTION 


1.1  THREE-DIMENSIONAL  EYE  ROTATIONS 

This  report  explains  the  basic  mathematics  necessary  to  understand  three-dimensional  (3D) 
eye  rotations,  in  the  easiest  and  the  simplest  way  it  can  be  presented. 

The  mathematical  background  assumed  is  not  more  demanding  than  the  college  freshman 
level. 

The  goal  of  this  report  is  to  show  easy-to-understand  derivations  for  the  key  equations  used 
in  the  mathematics  for  3D  eye  rotations.  This  is  based  on  the  philosophy  that  the  ideal  way  to  under¬ 
stand  and  apply  the  equations  is  to  be  able  to  understand  their  derivations.  To  the  author’s  best 
knowledge  at  this  writing,  there  seems  to  be  no  report  or  book  of  this  nature  available.  This  writing 
is  a  modest  endeavor  toward  that  direction,  hoping  someday  someone  will  take  over  and  continue 
the  task. 

Analysis  of  3D  eye  rotations  requires  a  lot  of  mathematics,  and  the  particular  mathematical 
methods  that  are  commonly  used  stem  from  many  areas  that  are  themselves  a  wide  field.  This  report 
only  touches  on  the  standard  methods  used  today.  Many  other  methods  are  used  and  can  be  used 
to  gear  the  analysis  to  larger  movements  (such  as  the  head  movements).  The  reader  must  keep  in 
mind  that  every  method  has  its  limitations  for  eye  movement  analysis  and  should  be  used  cautiously. 

Section  1  introduces  virtually  all  the  basic  linear  algebra,  along  with  Appendix  A,  required  to 
understand  this  report.  Section  2  describes  and  demonstrates  Euler’s  Theorem,  which  is  a  central 
concept  in  studying  3D  rotations.  Section  3  describes  a  geometrical  approach  to  understanding 
Listing’s  and  Donders’  Laws.  These  laws  are  immediate  applications  of  this  theorem  in  eye 
movement  analyses. 

Section  4  covers  Euler  Angles.  These  are  the  minimum  numbers  of  parameters  required  to 
describe  a  rotation  in  space.  Two  sets  of  Euler  Angles  used  most  commonly  in  eye  movement 
applications,  namely  Fick’s  and  Helmholtz’s  coordinate  systems,  are  described  in  this  section. 

Central  to  understanding  an  application  of  the  3D  rotations  is  the  concept  of  Rotation 
Matrices  which  is  covered  in  Section  4.  For  a  person  dealing  with  3D  rotations,  it  is  important  to 
understand  rotation  matrices,  how  to  generate  them  within  an  experimental  setup,  and  what  they 
represent.  Sections  5  through  13  include  the  mathematical  and  geometrical  descriptions  of  these 
matrices,  their  characteristics,  and  their  relationships  to  Fick’s  and  Helmholtz’s  coordinates. 

Sections  14  and  15  establish  the  concept  of  angular  velocity  in  3D  rotations  by  demonstrating 
the  Theorem  of  Coriolis. 
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Section  16  and  the  remaining  sections  turn  our  focus  to  another  method  of  describing  3D 
rotations.  As  a  result  of  Euler’s  Theorem,  instead  of  describing  a  rotation  in  space  as  three  consecu¬ 
tive  rotations  (resulting  in  Euler  angles)  we  can  define  a  set  of  parameters  which  defines  an  equiva¬ 
lent  single-axis  of  rotation  and  the  angle  of  rotation  around  this  axis.  Quaternions  and  rotation 
vectors  are  the  two  most  widely  used  of  such  methods  for  describing  rotations.  These  methods  and 
their  special  mathematics  are  described  in  Sections  16  through  26. 


1.2  REVIEW  OF  SOME  ALGEBRA 

To  start,  we  review  the  Pythagorean  Theorem.  Referring  to  Figure  1-1,  we  see  that 
cos  0  =  x/r  where  r  is  the  magnitude  of  the  vector  r,  or 

x  =  rcos0  (1-1) 

which  is  the  component  of  r  along  the  X-axis. 

Also,  we  see  that  sin0  =  y/r,  or 

y  =  rsin0  =  rcos  [90°  —  0]  (1-2) 

which  is  the  component  of  r  along  the  Y-axis. 


The  method  of  determining  the  components  of  a  vector  (such  as  r  in  Figure  1-1)  along  the 
coordinate  axes  this  way  is  identical  to  the  method  for  determining  the  elements  of  a  Rotation 
Matrix,  or  Coordinate  Transformation  Matrix  (see  Section  5,  the  Rotation  Matrix). 


Y 

i 

i 

v 

P  (x,y) 

J 

A 

J 

fc-  V 

0 

i  } 

C 

Figure  1-1.  Example  of  Pythagorean  Theorem 
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From  Equation  (1-1)  and  Equation  (1-2): 

x2  +  y2  =  [rcos0]2  +  [rsin0]2 
=  r2  [cos20  +  sin2  0] 

=  r2  (since  cos20  +  sin20  =  1)  (1-3) 

which  is  the  familiar  Pythagorean  Theorem,  and  validates  indirectly  the  method  used  to  determine 
the  x  and  y  components  of  r.  Denoting  the  unit  vector  along  the  X-axis  by  i,  and  along  the  Y-axis 
by  j,  we  can  represent  the  vector  r  in  terms  of  x  and  y  by: 

r  =  ix  +  jy  or  r  =  jy  +  lx  (1-4) 

Note  in  Equation  (1-4)  that  the  vector  addition  commutes.  We  can  do  this  because  (referring 
to  Figure  1-1)  we  can  reach  the  point  p  (x,  y)  either  by  moving  to  the  right  first  and  then  upward, 
or  by  moving  upward  first  and  then  to  the  right.  We  emphasize  this  because  the  consecutive  angular 
motions  generally  do  not  commute.  This  will  be  demonstrated  later. 


1.3  DOT  PRODUCT 


Denote  two  vectors  A  and  B  in  a  3D  orthogonal  frame  by 
A  =  axi  +  ayj  +  azk 

B  =  bxi  +  byj  +  bzk 


where  i,  j,  and  k  are  the  unit  vectors  along  the  x,  y,  and  z  axes  respectively. 

The  Dot  Product  (also  called  Scalar  Product  or  Inner  Product)  denoted  by  A  •  B  (pro¬ 
nounced  as  A  dot  B)  is  defined  to  be 

A  •  B  =  |A||B|cos0  (1-7) 

where  0  is  the  angle  between  A  and  B.  Noting  that  |i|  =  |j|  =  1  and  ilj,  we  have  using  Equa¬ 
tion  (1-7): 


i  •  j  =  cos  |  =  0 


(1-8) 
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Similarly, 


ij=j  k  =  k-  i  =  cos  ^  =  0 
j  i  =  k  j  =  i  k  =  cosf  =  0 

also, 

i  •  i  =  j  •  j  =  k  •  k  =  cos  0  =  1 
Therefore, 

A  •  B  =  [axi  +  ayj  +  azk]  •  [bxi 


(1-9) 

d-10) 

(1-11) 

byj  +  bzk]  (1-12) 


reduces  to,  using  Equations  (1-9),  (1-10),  and  (1-11): 

A  ■  B  =  axbx  +  ayby  -I-  azbz  =  B  •  A 


(M3) 


The  result  of  Equation  (1-13)  is  a  scalar  value  describing  the  projection  of  vector  A  onto 
vector  B  (or  vice  versa).  If  A  and  B  are  orthogonal  to  each  other  (0=90°),  then  A  •  B  =  0. 

Expressing  A  and  B  by  column  vectors, 


A  = 

V 

ay 

and  B  = 

n>*l 

by 

bz 

we  may  express  A  •  B  equivalently  by 


AtB  =  [ax  ay  az] 


-  axbx  +  ayby  +  azbz  (1-15) 

The  transpose  of  A,  denoted  by  AT ,  is  obtained  by  converting  the  columns  of  A  into  the  rows 
of  A  one  at  a  time  in  sequence. 

Equation  (1-15)  above  yields  the  same  results  as  A  •  B  given  in  Equation  (1-13).  ATB  is 
referred  to  as  the  Inner  Product,  while  A  •  B  is  referred  to  as  the  Dot  Product  although  they  both 
mean  the  same  thing.  It  is  also  called  Scalar  Product  because  the  results  are  scalar,  not  vector. 
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1.4  CROSS  PRODUCT 

The  Cross  Product  (also  called  Vector  Product)  of  two  vectors  A  and  B  denoted  by  A  x  B 
(pronounced  as  A  cross  B)  is  another  vector,  defined  by: 


A  x  B  =  |A||B|sin0u 


d-16) 


where  0  is  the  angle  between  vectors  A  and  B,  and  u  is  a  unit  vector  perpendicular  to  the  plane  con¬ 
taining  both  A  and  B,  and  the  direction  of  u  is  along  the  right  hand  thumb  if  its  fingers  curl  from 
A  to  B.  Hence,  this  convention  is  called  right  hand  rule. 

Note  the  result  ofEquation  (1-16)  is  vector  orthogonal  toboth  vectors  AandB.  If  A  x  B  =  0, 
then  it  means  the  vectors  A  and  B  lie  along  the  same  line  implying  0  must  be  0. 

It  follows,  using  Equation  (1-16): 

i  x  j  =  sin^k  =  k;  j  X  i  =  -k 

j  x  k  =  sin^i  =  i;  k  x  j  =  -i 


k  x  i  =  sin|j  =  j;  i  x  k  =  -j 


(1-17) 


Note 


i  x  i  =  sin(0)k  =  0  [zero  vector] . 


Similarly, 

jxj=kxk  =  ixi  =  0  [zero  vector] . 

It  follows  using  Equation  (1-5)  and  Equation  (1-6): 

A  x  B  =  [axi  +  ayj  +  azk]  x  [bxi  +  byj  +  bzk] 

which  reduces  to,  using  Equations  (1-16)  to  (1-19): 

A  x  B  —  (aybz  azby)i  +  (azbx  axbz)j  +  (axby  aybx)k 

Now,  the  following  determinant: 

|i  j  k 

ax  ay  az 

|bx  by  bz. 


(1-18) 


(1-19) 


(1-20) 


(1-21) 


(1-22) 
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gives  the  same  results  as  Equation  (1-21),  and  is  much  easier  to  compute.  Therefore,  A  x  Bis  con¬ 
ventionally  expressed  by: 


A  x  B 


i  j  k 

ax  %  az 

bx  by  bz 


(1-23) 


and  can  be  treated  as  the  definition  of  the  Cross  Product  given  in  Equation  (1-16). 


1.5  MATRIX  REPRESENTATION  OF  THE  CROSS  PRODUCT 


For  R  =  rxi  +  ryj  +  rzk  and  W  =  wxi  +  wyj  +  wzk,  we  get,  using  Equation  (1-21)  and 
Equation  (1-23): 


W  x  R 


i  j  k 

WX  Wy  WZ 

fx  ry  rz 


=  [wyrz  -  wzry]i 


+  [wzrx  -  wxrz]j 


+  [wxry  -  wyrx]k 


(1-24) 


Now  represent  the  vector  R  =  rxi  +  ryj  +  rzk  by  an  equivalent  Column  Matrix 


and  a  vector  W  =  wxi  +  wyj  +  wzk  by 


W  = 


wx 

Wy 

Wz 


We  want  to  find  a  matrix  representation  W*R  corresponding  to  the  vector  representation 
x  R  given  in  Equation  (1-24)  such  that  the  components  of  W*R  is  equal  to  the  components  of 
x  R. 


1-6 


NSWCDD/MP-99/17 


It  turns  out  that  if  we  express  W*  by: 


W*  = 


0 

w7 


-w. 


-w2 

0 

Wx 


Wv 

-WX 

0 


(1-25) 


and  perform  the  operation  W  R  ,  we  have: 

Wv 


W  R  = 


■  0 

— wz 

wz 

0  - 

& 

i  , 

Wx 

"wy  rz 

-  wzry" 

wzrx 

-  wxrz 

wxry 

X 

u 

£ 

1 

0 


(1-26) 


which  is,  component-wise,  equal  to  Equation  (1-24). 


It  follows  that  the  vector  Cross  Product  W  x  R  may  be  equivalently  expressed  by  the  Matrix 
Product  of  W*R  in  the  sense  that  its  components  are  the  same.  This  is  convenient  because,  while 
W  x  R  gives  better  physical  insight,  W*R  is  much  easier  in  computer  applications. 

If  we  take  the  transpose  of  Equation  (1-25): 


I W 


nT 


‘  0 

— wz 

Wy- 

= 

wz 

0 

-Wx 

-Wy 

Wx 

0 

0 

Wz 

-Wy" 

= 

-w2 

0 

Wx 

_Wy 

-Wx 

0 

'  0 

— wz 

Wy  ‘ 

zrz  _ 

wz 

0 

-Wx 

-Wy 

wx 

0 

=  -W 

We  find  that  matrix  W*  is  skew-symmetric  because 
[W*]T  =  -W*  or  W*  =  -[W*]T. 


(1-27) 


(1-28) 
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See  the  following  paragraph. 

1.6  SOME  DEFINITIONS  OF  MATRICES 

The  matrix  obtained  by  interchanging  the  rows  and  columns  of  matrix  A  one  by  one  is  called 
the  Transpose  of  A,  and  is  denoted  by  AT  (A  transpose). 


For  example,  for 
A  = 


1  2  3 

4  5  6 


,  A 


T  _ 


'1  4' 

2  5 

3  6 


[Af  - 


'1  4' 

T 

2  5 

3  6 

= 

’l  2  3’ 
4  5  6 

It  follows  that. 


[Af  =  A 


(1-29) 


It  can  be  easily  shown  by  direct  substitutions  (of,  say,  three-by-three  matrices)  that, 

(A  +  B)t  =  At  +  Bt  (1-30) 

(kA)T  =  kAT  (where  k  is  a  scalar.)  (1-31) 

(AB)T  =  BTAT  (see  Appendix  A,  Equation  (A-8))  (1-32) 


A  Square  Matrix  such  that, 

A  =  AT 

is  called  a  Symmetric  Matrix.  Thus,  for  a  three-by-three  Symmetric  Matrix  A: 


0-33) 


A  = 


lll 

a12 

a13 

l21 

a22 

a23 

l31 

a32 

a33 

lll 

a21 

a31 

l12 

a22 

a32 

l13 

a23 

a33 

=  At 


d-34) 
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A  Square  Matrix  B  is  the  inverse  of  a  Square  Matrix  A  if 
AB  =  BA  =  I 


d-35) 


The  inverse  of  A,  if  it  exists,  is  denoted  by  A  Thus  using  B=A  1  in  Equation  (1-35), 

AA"1  =  A-1  A  =  I  (1-36) 

where  I  is  the  Identity  Matrix.  For  instance,  the  three-by-three  Identity  Matrix  is  given  by: 


I  = 


1  0  O’ 
0  1  0 
0  0  1 


(1-37) 


A  Square  Matrix  A  is  called  an  Orthogonal  Matrix  if 
AAt  =  AtA  =  I 

Thus,  for  a  three-by-three  Orthogonal  Matrix  A: 

AAt  = 


(1-38) 


'all 

a12 

a13* 

‘an 

a21 

a3l" 

a21 

a22 

a23 

a12 

a22 

a32 

a31 

a32 

a33 

a13 

a23 

a33 

‘all 

a21 

2 

l3l" 

'all 

a12 

a13* 

a12 

a22 

a32 

a21 

a22 

a23 

a13 

a23 

2 

l33 

a31 

a32 

a33 

T 

0 

O' 

ata 

=  I 

= 

0 

1 

0 

0 

0 

1 

(1-39) 

Comparing  Equation  (1-36)  with  Equation  (1-38),  we  conclude  the  following  for  an 
Orthogonal  Matrix  A, 

At  =  A"1  (1-40) 

It  is  demonstrated,  in  Appendix  A,  using  two  three-by-three  Orthogonal  Matrices  A  and  B 

that, 

[AB]T  =  BtAt  =  B-1  A"1  =  [AB]-1  (1-41) 


l-9/(l-10  blank) 


NSWCDD/MP-99/17 


SECTION  2 

EULER’S  THEOREM  OF  A  SINGLE,  EQUIVALENT  ROTATION 


Euler’s  Theorem  states  that  if  a  moving  frame  initially  coincident  with  a  fixed  reference  frame 
makes  any  number  of  rotations,  regardless  of  how  it  reaches  the  final  orientation,  there  always  exists 
a  single  equivalent  rotation  with  a  finite  angle  about  a  single-axis  through  the  origin. 

In  Section  3,  we  will  describe  a  direct  application  of  this  theorem  to  eye  rotation  analysis. 

To  demonstrate  the  theorem,  consider  an  orthogonal  Frame  B  (with  axes  Xg,  Yb,  and  Zg)  rotat¬ 
ing  about  an  arbitrary-axis  u  fixed  in  space  (called  the  rotation-axis)  relative  to  another  orthogonal 
Frame  A  (with  axes  Xa,  Ya  and  Za),  which  is  fixed  in  space  and  serves  as  a  reference  frame. 
Frame  B  is  initially  coincident  with  Frame  A,  as  shown  in  Figure  2-1. 


Zb 

zA 

ROTATION-AXIS 

— -S^ZA’  0ZB  u 

0Ya,0YB 

1  ya,  Yg 

®XA®XB 

xA>  xB 

INITIALLY  FRAME  B  IS  COINCIDENT  WITlf  FRAME  A 

Figure  2-1.  Example  of  Euler’s  Theorem 


In  Figure  2-1,  u  is  an  unit  vector  rigidly  attached  to  the  movable  Frame  B,  while  maintaining 
a  fixed  direction  relative  to  Frame  A,  (but  not  attached  to  Frame  A).  The  unit  vector  u  makes 
angle  0x  with  both  Xa  and  Xg  axes,  angle  0y  with  both  Ya  and  Yg  axes,  and  angle  Qz  with  both 
ZA  and  Zg  axes  initially  as  shown  in  Figure  2- 1 .  Note  that  the  angles  0xa,  ©ya.  and  0za  that  u  makes 
with  Frame  A  remain  constant,  because  the  direction  of  u  relative  to  Frame  A  is  fixed.  Also,  the 


2-1 


NSWCDD/MP-99/17 


angles  0xb,  9yb,  and  Qzb  that  u  makes  with  Frame  B  remain  constant,  because  Frame  B  is  rotating 
around  the  axis  of  u,  and  u  is  rigidly  attached  to  Frame  B. 


Now,  rotate  Frame  B  around  u  by  an  arbitrary  angle  4>  while  maintaining  the  direction  of  u 
fixed  relative  to  Frame  A  (relative  to  space),  as  shown  in  Figure  2-2. 


Figure  2-2.  Frame  B  Rotated  Around  u 


The  vector  u  is  located  in  the  reference  Frame  A  by  the  angles  0XA,  0YA,  and  0ZA,  as  shown 
in  Figure  2-2.  After  the  rotation,  u  is  located  in  the  rotating  Frame  B  by  the  angles  0XB,  Oyg,  and 
0ZB.  Because  of  Euler’s  theorem,  and  in  reality, 

®XA  =  ®XB,  ®YA  =  ®YB,  ®ZA  =  ®ZB  (2-1) 
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We  reached  the  orientation  of  Frame  B  in  Figure  2-2  by  an  arbitrary  rotation  <j)  around  a  unit 
vector  pointing  to  an  arbitrary  (but  fixed)  direction.  This  rotation  does  not  change  0XB,  Qyg,  or 
0ZB,  as  previously  explained. 

The  Frame  B  that  was  initially  aligned  with  Frame  A  could  have  reached  the  final  orientation 
of  Frame  B  in  Figure  2-2  by  a  sequence  of  many  rotations  about  many  different  axes  of  rotation. 
But  we  are  searching  for  a  unique  way  to  define  the  orientation  of  Frame  B  relative  to  Frame  A. 
Establishing  u  by  0X,  0Y,  and  0Z  and  rotation  angle  <f)  produces  this  unique  definition  of  the 
orientation  of  Frame  B  relative  to  Frame  A. 

Therefore,  it  is  obvious  that  for  any  orientation  of  Frame  B ,  there  always  exists  a  single,  equiva¬ 
lent  rotation-axis  with  an  equivalent  rotation  angle  that  could  move  Frame  B  from  its  initial  orienta¬ 
tion  in  the  reference  Frame  A  to  its  final  orientation.  As  an  analogy,  consider  two  points  in  a  plane, 
such  as  P  and  Q.  Since  there  are  many  paths  that  go  from  P  to  Q,  we  want  to  define  a  single  path 
—  the  straight  line  —  as  the  unique  passage  from  P  to  Q. 

Remember,  unit  vector  u  always  makes  the  same  orientation  angles  0x,  By,  and  0z  with  respect 
to  both  Frame  A  and  Frame  B.  Thus,  it  has  the  same  components  in  both  Frame  A  and  Frame  B. 
This  is  shown  in  the  following,  denoting  the  magnitude  of  u  by  u  which  is  unity: 

xA  =  ucos0x  =  cos0x  =  xB 
yA  =  u  cos  0y  =  cos0y  =  yB 

zA  =  ucos0z  =  cos  0z  =  zB  (2-2) 
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SECTION  3 

LISTING’S  LAW  AND  DONDERS’  LAW  BY  SIMPLE  GEOMETRY 

Euler’s  Theorem  of  single  equivalent  rotation  has  an  immediate  application  to  eye  movements 
expressed  in  Donders’  Law  and  Listing’s  Law. 

Donders’  Law  states  that  every  eye  position  in  space  is  described  by  a  single  orientation  (fixed 
bearing  and  elevation  -  a  and  (3  in  Figure  3-1)  of  the  eye. 

Listing’s  Law  states  that  the  eye  can  reach  the  final  orientation  through  a  single  equivalent 
rotation  about  an  axis  of  rotation.  The  axes  of  rotation  for  all  final  eye  orientation  will  be  in  the  same 
plane,  called  Listing’s  Plane. 


P  FINAL  POSITION 
>  (LINE  OF  SIGHT) 


Y  (LEFT) 


xt 


X  (FORWARD) 


INITIAL  POSITION 
(LINE  OF  SIGHT) 

PRIMARY/REFERENCE 

POSITION 


Figure  3-1.  The  Orientation  of  Line  of  Sight 


Referring  to  Figure  3-1,  the  line  of  sight  of  the  initial  eye  position  is  OX  along  the  X-axis, 
which  we  may  consider  as  the  primary/reference  position;  its  final  position  is  OP  with  the  bearing 
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a  and  the  elevation  (3.  The  eye  may  have  reached  the  final  position  OP  from  the  primary  position 
by  various  sequences  or  paths. 

Recall  from  Euler’s  Theorem  that  there  exists  an  axis  of  rotation  which  brings  OX  to  OP  by 
a  single  equivalent  rotation,  regardless  of  how  the  final  position  OP  may  have  been  reached  from 
the  initial  position  OX.  Listing’s  Law  uses  Euler’s  Theorem  to  find  the  Listing  plane. 

We  now  describe  two  methods  for  finding  Euler’s  axis  of  single  equivalent  rotation. 

In  the  first  method,  Euler’s  axis  of  rotation  is  found  by  constructing  a  perpendicular  line, 
through  origin,  to  the  plane  defined  by  lines  formed  by  OX  and  OP.  The  perpendicular  line  is  OZ' 
in  Figure  3-2  and  OZ"  in  Figure  3-3.  The  line  is  the  Euler’s  rotation-axis  that  brings  the  OX  to  the 
OP  by  the  right  hand  rule,  with  the  thumb  pointing  in  the  direction  of  the  perpendicular  line,  OZ'. 

Denoting  the  unit  vector  along  OX  by  ux  and  the  unit  vector  along  OP  by  up ,  the  vector  cross 
product  ux  X  Up  yields  a  vector  L  which  is  perpendicular  to  the  plane  formed  by  ux  and  up.  That 
is. 


L  =  ux  x  up  (3-1) 

Referring  to  Figure  3-2,  expressing  up  by  its  components  in  the  X-Y-Z  frame  by 

Up  =  ixp  +  jyp  +  kzp  (3-2) 


and  realizing  that  ux  is  equal  to  i  since  it  is  a  unit  vector  along  the  X-axis,  we  have  based  on 
Equation  (3-1),  using  Equations  (1-17)  and  (1-18): 

L  =  i  x  (ixp  +  jyp  +  kzp) 

=  i  x  ixp  +  i  x  jyp  +  i  x  kzp 
=  0  +  kyp  -  jzp 

=  jzP  +  kyP  (3-3) 


Equation  (3-3)  shows  that  L  must  lie  in  the  Y-Z  plane  since  it  has  no  i  component.  Since  up 
may  represent  any  final  eye  position,  the  Equation  (3-3)  shows  that  all  axes  of  rotation  corresponding 
to  any  final  eye  position  must  lie  in  the  Y-Z  plane,  which  is  called  Listing’s  Plane. 

The  Equation  (3-3)  suggests  that  the  angle  <f>  which  L  makes  with  the  Z-axis  (see  Figure  3-2) 
may  be  computed  from: 


tancj) 


ypi 


or  <j)  =  tan 


(3-4) 
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Referring  to  Figure  3-1,  if  we  know  the  bearing  a  and  the  elevation  (3  of  OP  or  up ,  we  may 
determine  xp,  yp,  and  zp  in  Equation  (3-2)  by  inspection  as: 

xp  =  (cos  P)  cos  a 
yp  =  (cos  P)  sin  a 

zp  =  sin  p  (3-5) 

In  the  second  method,  Euler’s  rotation-axis  is  found  by  rotating  coordinate  frames  in  the  fol¬ 
lowing  manner  which  shows  that  the  angle  <j)  is  the  torsion  angle  about  the  X-axis  (primary/refer¬ 
ence  axis). 

Step-1;  Construct  a  plane  containing  OX  and  OP;  find  the  intersection  of  this  plane  with  the 
original  Y-Z  plane  as  shown  in  Figure  3-2.  Call  this  intersection  the  Y'-axis.  The  line  OY'  makes 
an  angle  <|)  with  the  OY-axis,  or  ZYOY'=<|). 


Z 


Figure  3-2.  Rotation  of  X-Y-Z  Frame  about  X-Axis  (Primary/Reference  Axis) 

Step  2:  Rotate  the  X-Y-Z  frame  about  the  X-axis  until  the  Y-axis  reaches  the  Y'-axis,  or  until 
the  rotation  angle  about  X-axis  is  equal  to  <j).  Call  the  new  frame  X'-Y'-Z'. 


3-3 


NSWCDD/MP-99/17 


Then  X-axis  coincides  with  X'-axis,  while  the  Y-axis  moves  to  Y'-axis  through  angle  (|).  The 
Z-axis  moves  to  Z -axis  also  with  angle  <J)  while  staying  on  the  original  Y-Z  plane. 

Note  that  both  Z-axis  and  Y'-axis  are  perpendicular  to  the  X-axis  because  they  remain  in  the 
original  Y-Z  plane. 

Step  3:  Now  rotate  X'-Y'-Z'  frame  about  the  Z'-axis  by  the  angle  0,  until  OX  (primary  posi¬ 
tion)  coincides  with  OP  (final  position)  as  shown  in  Figure  3-3.  Call  this  new  frame  X"-Y"-Z" 
frame. 


Then  the  X-axis  (or  the  X'-axis)  moves  to  X"-axis,  Y'-axis  moves  to  Y"-axis,  and  Z'-axis 
remains  at  the  same  direction  and  becomes  Z"-axis. 
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Note  that  the  Z"-axis  (Z-axis)  lies  on  the  original  Y-Z  plane.  The  X"-axis  lies  on  the  X'-Y' 
plane  (or  X'-Y'  plane).  The  Y"-axis  is  located  behind  the  vertical  Y-Z  plane,  while  still  lying  on 
the  X'-Y'  plane.  The  angle  0  is  the  single,  equivalent  rotation  angle  that  makes  the  primary  eye  posi¬ 
tion  OX  to  the  final  eye  position  OP  or  OX". 

To  summarize,  in  order  to  move  the  eye  from  the  primary  eye  position  OX  (X-axis)  to  the  final 
eye  position  OP  (X"-axis)  by  a  single  equivalent  rotation,  the  eye  has  to  be  rotated  about  the  Z"-axis 
in  Listing’s  plane  (which  is  the  Y-Z  plane)  by  angle  0  (the  single  equivalent  rotation).  Regardless 
of  any  direction  of  the  line  OP,  we  can  see  that  the  axis  of  single  equivalent  rotation  always  lies  on 
the  unique  plane  called  the  Listing’s  plane. 

We  see  that  Donders’  Law  and  Listing’s  Law  imply  that  the  angular  rotation  (J)  (in  Figures  3-2 
and  3-3),  called  cyclotorsion  or  simply  torsion,  is  fixed  for  each  final  position  (line  of  sight),  regard¬ 
less  of  the  different  paths  the  eye  might  have  taken  to  reach  there. 

So,  practically,  we  can  test  Listing  Law  by  measuring  torsion  <j)  at  any  final  eye  position.  The 
dot  product  of  the  corresponding  axis  of  rotation  with  the  primary/reference  position  should  be  equal 
to  0  (perpendicular). 

In  reality,  the  primary  position  is  slightly  lifted  upward  from  the  head’s  straight  forward  refer¬ 
ence  position.  This  is  equivalent  to  the  rotation  of  the  initial  eye  position  (which  is  coincident  with 
the  reference  head  position)  by  a  small  angle  around  the  Y-axis. 

For  this  situation,  the  Listing’s  plane  is  tilted  backward  by  the  same  small  angle,  and  its  axis 
may  be  determined  in  the  displaced  eye  frame  by  the  similar  method  as  explained  in  this  section. 
Coordinates  may  be  transformed  back  to  the  reference  head  frame  by  means  of  a  Rotation  Matrix 
for  the  Y-axis  rotation.  This  will  be  discussed  in  Section  5. 
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SECTION  4 
EULER  ANGLES 

Euler  Angles  orient  one  orthogonal  Frame  B  displaced  relative  to  another  orthogonal  Frame  A 
by  making  three  sequential  rotations  of  Frame  B  relative  to  Frame  A.  The  Frame  B  is  initially 
aligned  with  the  Frame  A  by  a  common  origin.  The  first  rotation  is  about  any  axis  of  the  initial 
Frame  B.  The  second  rotation  is  about  either  of  the  two  axes  of  the  displaced  Frame  B  not  used  for 
the  first  rotation.  The  third  rotation  is  about  either  of  the  two  axes  of  the  Frame  B  not  used  for  the 
second  rotation.  Thus,  the  number  of  permutations  for  the  possible  sequences  of  three  Euler  Angles 
is  (3)  (2)  (2)  or  12. 

One  particular  sequence  or  set  of  Euler  Angles  has  had  wide  application.  In  this  set,  the  initial 
Frame  Bo  with  Xo,  Yo,  and  Zo  axes  is  rotated  about  its  Zo-axis  through  an  angle  ({),  resulting  in  the 
Frame  Bi  with  Xj,  Yi  and  T\  axes  (see  Figure  4-1).  Note  that  the  Zi-axis  coincides  with  the 
Zo-axis  during  the  first  rotation.  In  the  second  rotation,  the  Frame  Bi  is  rotated  about  its  Xi-axis 
through  an  angle  0,  resulting  in  the  Frame  B2  with  X  2,  Y2,  and  Z2  axes  (see  Figure  4-2).  Note  that 
the  X2~axis  coincides  with  the  Xi-axis  during  the  second  rotation.  In  the  third  rotation,  the 
Frame  B2  is  rotated  about  its  Z2-axis  (displaced  Zi-axis)  by  an  angle  ip,  resulting  in  the  Frame  B3 
with  X3,  Y3,  and  Z3  axes  (see  Figure  4-3).  Note  that  the  Z3-axis  coincides  with  the  Z2-axis  during 
the  third  rotation.  This  Zo,  Xi,  Z2  Euler  Set  is  commonly  used  for  robotic  applications,  and  for 
aircraft  and  missile  applications. 


Out  of  twelve  possible  sets  of  Euler  Angles,  only  two  sets  are  historically  used  in  eye  rotation 
analyses:  Fick’s  system  and  Helmholtz’s  system  of  Euler  Angles.  Fick’s  system  will  be  described 
in  Section  4.1,  and  Helmholtz’s  system  will  be  described  in  Section  4.2. 
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Figure  4-1.  Frame  Bi 
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Figure  4-3.  Frame  B3 
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4.1  FICK’S  SYSTEM  OF  EULER  ANGLES 

In  the  Fick’s  system,  as  in  the  Helmholtz’s  system,  initially  both  the  Head  frame  and  the  Eye 
frame  are  coincidentally  aligned  with  both  of  their  X-axes  pointed  straight  forward,  Y-axes  hori¬ 
zontally  leftward,  and  Z-axes  vertically  upward. 

The  sequence  of  rotations  in  Fick’s  system  is  about  the  Zo-axis  of  the  initial  Frame  Eo,  followed 
by  rotation  about  the  Y i-axis  of  the  Frame  E]  (the  displaced  Frame  Eq)  followed  by  rotation  about 
the  X2-axis  of  the  Frame  E2  (displaced  Frame  Ei).  This  is  a  Zo,  Yi,  X2  Euler  Set. 

For  the  first  rotation  of  Fick’s  system,  the  initial  Frame  Eo  (E  for  Eye)  with  its  Xq,  Yo,  Zo  axes 
is  rotated  about  the  Zo-axis  resulting  in  the  Frame  Ei  with  X\,Y\,Zi  axes.  This  is  shown  in 
Figure  4-4  using  the  rotation  angle  of  90°  for  clarity  and  ease  of  comprehension.  We,  of  course, 
understand  that  90°  eye  rotation  is  not  realizable  physiologically.  All  eye  rotation  angles  of  Fick’s 
system,  as  well  as  those  of  Helmholtz’s  system,  are  constrained  to  be  well  less  than  90°  in  reality. 

For  the  second  rotation,  the  Frame  Ei  is  rotated  about  its  Y i-axis,  resulting  in  the  Frame  E2 
with  X2,  Y2,  Z2  axes.  This  is  shown  in  Figure  4-5  using  the  rotation  angle  of  90°.  Note  that  Y i-axis 
is  pointed  into  the  paper.  Thus,  the  thumb  is  pointed  into  the  paper  as  well,  when  the  right  hand  rule 
is  applied  (Figure  4-5). 

For  the  third  rotation,  the  Frame  E2  is  rotated  about  its  X2-axis,  resulting  in  the  Frame  E3  with 
X3,  Y3,  Z3  axes  (Figure  4-6). 
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(B)  AFTER  THE  ROTATION 


POINTS  OUT  OF  THE  PAPER  LIKE  AN 
ARROW  HEAD 

POINTS  INTO  THE  PAPER  LIKE  AN 
ARROW  TAIL 


Figure  4-4.  First  Rotation  About  the  Zo-Axis  in  Fick’s  System 
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Figure  4-5.  Second  Rotation  About  the  Yj-Axis  in  Fick’s  System 
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Figure  4-6.  Third  Rotation  About  the  X2-Axis  in  Fick’s  System 
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4.2  HELMHOLTZ’S  SYSTEM  OF  EULER  ANGLES 

In  the  Helmholtz’s  system,  as  in  the  Fick’s  system,  initially  the  Head  frame  coincides  with  the 
Eye  frame  with  the  X-axis  pointed  forward,  the  Y-axis  leftward,  and  Z-axis  upward. 

The  sequence  of  rotations  in  Helmholtz’s  system  is  first  about  the  Yo-axis  of  the  initial 
Frame  Eo ,  followed  by  rotation  about  the  Zi-axis  of  the  Frame  Ei  (the  displaced  Frame  Eo),  fol¬ 
lowed  by  rotation  about  the  X2-axis  of  the  Frame  E2  (the  displaced  Frame  Ei ).  This  is  a  Y0,Zi, 
X2  Euler  Set. 

For  the  first  rotation  of  Helmholtz’s  system,  the  initial  Frame  Eq  is  rotated  about  its  Yo-axis, 
resulting  in  the  Frame  Ei  with  Xi,  Yi,  Z\  axes.  This  is  shown  in  Figure  4-7,  using  rotation  angle 
of  90°  again  for  demonstration  only.  All  rotation  angles  of  Helmholtz’s  system  are  constrained  to 
be  less  than  90°  in  reality. 

For  the  second  rotations,  the  Frame  Ei  is  rotated  about  its  Zi-axis  resulting  in  the  Frame  E2 
with  X2,  Y2,  Z2  axes.  This  is  shown  in  Figure  4-8,  using  a  rotation  angle  of  90°. 

For  the  third  rotation,  the  Frame  E2  is  rotated  about  its  X2-axis,  resulting  in  the  Frame  E3  with 
X3,  Y3,  Z3  axes.  This  shown  in  Figure  4-9,  using  the  rotation  angle  of  90°- 

Note  that  the  first  two  rotations  of  Fick’s  system  are  the  Zo-axis  rotation  followed  by  the 
Yi-axis  rotation,  while  those  of  Helmholtz’s  system  are  the  Yo-axis  rotation  followed  by  the 
Zi-axis  rotation,  which  is  the  reverse  order  of  the  former.  In  both  systems,  the  third  rotation  is  about 
the  X2-axis. 
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Figure  4-8.  Second  Rotation  About  the  Zi-Axis  in  Helmholtz’s  System 
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(B)  AFTER  THE  ROTATION 


Figure  4-9.  Third  Rotation  About  the  X2-Axis  in  Helmholtz’s  System 
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4.3  ANGULAR  ROTATIONS  DO  NOT  COMMUTE 

As  briefly  discussed  in  Section  1.2,  adding  two  vectors,  Vi  and  V2,  we  found  out  that  Vi  +  V2 
=  V2  +  Vi,  as  shown  in  Figure  4-10: 


Figure  4-10.  Addition  of  Two  Vectors 


That  is,  vector  motion  (addition)  commutes.  In  contrast,  the  first  two  rotations  of  Fick’s  sys¬ 
tem  (Zo-axis  rotation  followed  by  Y i-axis  rotation)  do  not  result  in  the  same  orientation  produced 
by  the  first  two  rotations  of  Helmholtz’s  system  (Yo-axis  rotation  followed  by  Zi-axis  rotation)  for 
the  same  amounts  of  the  angular  rotations,  as  repeated  below  in  Figure  4-11,  copied  from  Figures  4-5 
and  4-8. 
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Y2  =  Y, 


ORIENTATION  OF  EYE  FRAME  AFTER  THE 
FIRST  ROTATION  ABOUT  Z0-AXIS  ROTATION 
BY  90°  FOLLOWED  BY  THE  SECOND 
ROTATION  ABOUT  Y-AXIS  BY  90°  IN  FICK'S 
SYSTEM  (SEE  FIGURE  4-5). 


Y2 


ORIENTATION  OF  EYE  FRAME  AFTER  THE 
FIRST  ROTATION  ABOUT  Y0-AXIS  BY  90° 
FOLLOWED  BY  THE  SECOND  ROTATION 
ABOUT  Zf-AXIS  BY  90°  IN  HELMHOLTZ’S 
SYSTEM  (SEE  FIGURE  4-8). 


Figure  4-11.  Comparison  of  Eye  Frame  Orientation  of  Fick’s  System  and 
Helmholtz’s  System  After  the  Second  Rotation 


Mathematically,  this  means  that  for  any  Vectors  Vj  and  V2,  Vj  +  V2  =  V2  +  Vj.  But,  for 
Rotation  Matrices  Rj  and  R2,  generally  RjR2  ^  R2Ri-  That  is,  the  rotation  corresponding  to  Rj 
followed  by  the  rotation  corresponding  to  R2  is  not  the  same  as  R2  followed  by  Rj.  That  is,  Rj  and 
R2  do  not  commute.  This  is  discussed  in  Section  5. 

Although  finite  rotations  may  not  be  considered  as  vectors  and  thus  do  not  commute,  infinitesi¬ 
mal  rotations  may  be  considered  as  such.  That  is,  in  the  limit,  as  the  angle  of  rotation  becomes  very 
small,  Rotation  Matrices  commute  as  does  vector  addition.  This  is  discussed  in  Section  7. 
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SECTION  5 

THE  ROTATION  MATRIX 

In  the  following  pages,  we  will  derive  the  Rotation  Matrices  for  the  three  basic  rotations  around 
the  Z-axis  (yaw  or  horizontal  rotation),  the  Y-axis  (pitch  or  vertical  rotation),  and  the  X-axis  (roll 
or  torsional  rotation).  Any  3D  rotation  in  space  can  be  produced  by  a  combination  (multiplication) 
of  these  basic  rotations. 

5.1  BASIC  ROTATION  AROUND  THE  Z-AXIS 

Consider  a  vector  r  in  an  orthogonal  Frame  A  with  components  xA,  yA,  and  zA  as  shown  in 
Figure  5-1. 


Za 


Figure  5-1.  Frame  A  with  Vector  r 


Frame  B,  which  is  initially  coincident  with  Frame  A,  is  rotated  about  the  common  Za 
=  Zg-axis  counterclockwise  by  an  angle  0z  as  shown  in  Figure  5-2.  The  same  vector  r  has 
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components  xB,  yB,  and  zB  in  Frame  B.  Our  goal  is  to  find  a  relationship  between  xA,  yA,  za  and  xB, 
yB,  zB  of  the  same  vector  r. 

Let  rp  (in  Figure  5-4)  be  the  projection  of  the  vector  r  in  the  Y-Z  plane  of  both  frames. 
Decomposing  the  components  of  rp  in  Xa,  Ya,  and  Xb  directions,  we  get  Figures  5-3  and  5-4.  We 
have  a  pair  of  identical  right  triangles  with  hypotenuse  xa,  and  another  pair  with  hypotenuse  yA  in 
Figure  5-4. 

By  inspection,  we  get,  from  Figures  5-3  and  5-4: 

xB  =  xAcos0z  +  yAsin0z  (54) 

Since  yA  cos  0Z  =  yB  +  xA  sin  0Z,  we  get: 

Yb  =  Ya  cos  0z  ~  xa  sin  ez 

=  -  xA  sin  0Z  +  yA  cos  0Z  (5-2) 

Since  the  ZA-axis  coincides  with  the  ZB-axis, 


zB  =  zA  (5-3) 


Figure  5-2.  Frame  A  is  Rotated  to  Frame  B  about  the  Z-axis  by  Angle  0Z 
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Figure  5-3.  Components  of  r  in  Frame  A  decomposed 
into  components  in  Frame  B  (see  Figure  5-4) 
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Figure  5-4.  Components  of  rp  in  both  Frame  A  and  Frame  B 


5-4 


NSW  CDD/MP-99/17 


The  coordinates  Xa,  yA»  za  for  r  in  Frame  A  may  be  represented  either  as  a  column  vector 


yA 


or  as  row  vector  [xA  Ya  zb].  Similarly,  the  coordinates  xg,  yg,  zb  for  rin  Frame  B  may  be  represented 


either  as 


XB 

Yb 

ZB 


or  as  [xB  yB  zB]. 

If  we  use  column  vector  representation  for  components,  we  get  from  Equations  (5-1)  to  (5-3): 

cos  0Z  sin  0Z  0 

sin  0Z  cos  0Z  0  |  |  yA  |  (5-4) 

0  0  1 


XB 

Yb 

_ 

ZB 

XA 

Ya 

ZA 

By  denoting  column  vectors  by: 


rA  = 


XA 

Ya 

;  and 

rB  = 

XB 

Yb 

ZA 

ZB 

(5-5) 


and  denoting  the  Coefficient  Matrix  in  Equation  (5-4)  by  C?,  (from  Frame  A  to  Frame  B),  we  have 


cos  0Z 

sin  0Z  O’ 

Cj2 

^13 

— 

-sin  0Z 

cos  0Z  0 

C-21 

^22 

^23 

0 

0  1 

^31 

C32 

C33 

(5-6) 


It  follows: 


rB  = 


C?  rA 


(5-7) 


Superscript  B  in  rB  and  Superscript  A  in  rA  indicate  the  frame  in  which  we  have  resolved  the 
vector. 


The  matrix  CB  in  Equation  (5-7)  transforms  the  coordinate  of  a  vector  expressed  in  one  frame 
(Frame  A)  to  another  frame  (Frame  B).  Therefore,  it  is  called  the  “Coordinate  Transformation 
Matrix.”  It  is  also  called  “Rotation  Matrix”  because  it  shows  the  effect  of  rotation  of  one  frame 
(Frame  B)  relative  to  another  (Frame  A).  It  is  also  called  “Direction  Cosine  Matrix”  because  the 


elements  of  the  matrix  are  all  cosines  of  angle  0 


sin0  =  cos 


The  notation  CB  indi¬ 


cates  a  transformation  from  A  to  B,  or  from  the  subscript  to  the  superscript. 
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The  determinant  of  C®  in  Equation  (5-7)  denoted  by  |C®|  is 

I  cos0z  sin0z  0 

sin  0Z  cos  0Z  0 

0  0  1 

=  cos2  0Z  —  ( -  sin20z) 

=  1  (5-8) 

It  turns  out  that  the  determinant  of  any  Rotation  Matrix  is  always  equal  to  1. 

Reversing  the  rotation  process,  the  Frame  B  is  rotated  about  the  common  Za  =  Zg-axis  clock¬ 
wise  by  angle  0Z,  bringing  the  Frame  B  back  to  position  of  coincidence  with  Frame  A,  as  shown  in 
Figure  5-5,  thus  reversing  the  process  described  in  Figure  5-2. 


Figure  5-5.  Frame  B  Rotated  Clockwise  Back  to  Frame  A 


Knowing  0Z,  we  want  to  determine  xA,  yA,  zA  in  Frame  A  of  the  vector  r  in  terms  of  xB,  yB,  and 
zB  of  r  in  Frame  B. 
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Summing  components  in  Xa,  Ya  and  Za  directions: 


xA  =  XBCOS0Z  -  yBsin0z 

(5-9) 

yA  =  xB  sin  0Z  +  yBcos0z 

(5-10) 

ZA  =  ZB 

(5-11) 

Note  that  the  ( —  yB  sin  0Z)  term  in  Equation  (5-9)  is  negative  because  it  is  in  the  direction  of 
negative  XA-axis  direction. 

If  we  use  column  vector  representation  for  the  components,  we  get  from  Equations  (5-9), 
(5-10),  and  (5-11): 

XA  cos0z  -sin0z  0  xB 

yA  =  sin0z  cos  0Z  0  yB  (5-12) 

zA  0  0  1  zB 

Using  the  notations  given  in  Equation  (5-5),  and  denoting  the  Coefficient  Matrix  in  Equation  (5-12) 
by  CB,  (from  Frame  B  to  Frame  A),  we  have: 

cos0z  -sin0z  0 

C§  =  sin0z  cos  0Z  0  (5-13) 

0  0  1 

It  follows  that 

rA  =  C£  rB  (5-14) 

Substituting  Equation  (5-7)  into  Equation  (5-14)  for  rB,  we  have 

rA  =  c£  rA  (5-15) 

Since  rA  is  equal  to  itself,  Equation  (5-15)  indicates  that  CB  CB  must  be  the  Identity  Matrix: 

'1  0  O' 

CA  CB  =  I  =  0  1  0  (5-16) 

B  A  0  0  1 

This  makes  sense  geometrically  because  Frame  A  is  rotated  to  Frame  B  by  CB,  and  then 
Frame  B  is  rotated  back  to  Frame  A  by  CB.  Therefore,  Frame  B  is  identical  to  or  coincident  with 
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Frame  A.  Note  that  the  second  rotation  denoted  by  Cg  is  multiplied  from  the  left.  This  occurs 
because  the  coordinates  of  the  vectors  are  represented  in  column  vector  format  as  shown  in 
Equation  (5-5).  We  may  confirm  the  validity  of  Equation  (5-15)  analytically  by  direct  substitution: 


Using  Equation  (5-13)  and  Equation  (5-6): 


cos  0Z  -  sin  0Z  0 
sin  0Z  cos  0Z  0 
0  0  1 


cos0z  sin  0Z  0 
-sin  0Z  cos  0Z  0 
0  0  1 


cos20z  +  sin2  0Z  cos0z  sin0z  -  sin0z  cos0z 

sin  0Z  cos0z  -  cos0z  sin0z  sin2  0Z  +  cos2  0Z 
0  0 


0 

0 

1 


1  0  0 
0  1  0 
0  0  1 


(5-17) 


which  validates  the  conclusion  of  Equation  (5-16). 

Referring  C®  given  Equation  (5-6),  if  we  exchange  its  rows  and  its  columns,  that  is,  if  we  take 

T  T 

the  transpose  of  C®,  then  by  denoting  the  transpose  of  C®  by  (C®)  ,  we  find  (C®)  =  Cg.  That 
is, 

[c?]T  - 


sin  0Z  0* 
cos  0Z  0 
0  1 

=  C§  (5-18) 

which  is  the  same  as  Cg  given  in  Equation  (5-13). 

Substituting  Equation  (5-18)  into  Equation  (5-16), 


cos  0Z  - 
sin  0Z 
0 


cos0z  sin  0Z  0 

-sin0z  cos0z  0 
0  0  1 
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Now  substituting  Equation  (5-14)  into  Equation  (5-7): 

rB  =  CB  c£  rB  (5-20) 

Since  rB  is  equal  to  itself,  Equation  (5-20)  requires: 

c£  c£  =  I  (5-21) 

Now,  substituting  Equation  (5-18)  into  Equation  (5-21)  for  C£: 

Cl  [C®]T  =  I  (5-22) 

By  definition,  the  inverse  of  C  denoted  by,  C  “ 1  =  (CB )  1  must  satisfy  C  “ 1 C  =  CC  ~ 1  =  I, 
or  it  must  satisfy 

[CjT'cX-CSfcJf'-I  (5-23) 

From  Equation  (5-19)  and  Equation  (5-22): 

[cjf  Cl  =  Cl  [c>]T  =  I  (5-24) 

Comparing  Equation  (5-23)  and  Equation  (5-24): 

[ca]T  =  [ca]_1  (5-25) 

Equation  (5-25)  states  that  the  inverse  of  CB  is  equal  to  the  transpose  of  CB. 

For  any  Rotation  Matrix  C  that  transforms  the  coordinate  of  one  orthogonal  frame  to  another 
orthogonal  frame,  it  is  always  true  that 

CT  =  C"1  (5-26) 


A  matrix  that  satisfies  Equation  (5-26)  is  called  an  Orthogonal  Matrix. 

This  is  very  convenient  because  the  determination  of  C  -1,  which  is  generally  complicated,  may 
be  computed  simply  by  exchanging  the  rows  and  columns  of  C,  or  by  flipping  the  matrix  around  the 
diagonal-axis  for  an  Orthogonal  Matrix. 

From  Equation  (5-24)  we  have,  using  C  for  CB: 

CTC  =  I  =  CCT  (5-27) 
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Using  the  notations  defined  by  the  right  side  of  Equation  (5-6)  for  CB,  we  have 


C11 

^12 

C13 

C11 

^21 

C31 

CCT  = 

C21 

C22 

^23 

^12 

C22 

^32 

C31 

^32 

^33 

^13 

C23 

C33 

(5-28) 


Equation  (5-27)  or  Equation  (5-28)  is  called  the  orthogonality  condition.  Although  Equa¬ 
tion  (5-27)  is  derived  from  a  special  case,  it  is  true  for  all  Orthogonal  Matrices.  It  is  used  to  normal¬ 
ize  and  orthogonalize  the  C  matrix,  in  case  the  elements  Qj  determined  by  measurements  or  com¬ 
putations  deviate  from  the  condition  specified  by  Equation  (5-27)  or  Equation  (5-28).  In  that  case, 
this  is  geometrically  equivalent  to  non-orthogonal  frames.  For  example,  Equation  (5-28)  is  often 
used  to  orthogonalize  the  eye  coil  field  and  to  correct  for  measurement  errors.  (See  Sections  11  to 
13.) 


5.2  Z-AXIS  ROTATION  USING  ROW  VECTORS 


If  we  represent  the  coordinate  of  the  vector  (r A)T  by  a  row  vector  [xa  yA  za1  in  Frame  A  and 
the  coordinates  of  (rB)T  by  a  row  vector  [xb  yB  zb1  in  Frame  B,  we  get  from  Equations  (5-1)  to  (5-3): 


XB  Yb  zb]  =  [XA  yA  za] 


COS0Z 

sin0z 

0 


-sin0z  0 
cos  0Z  0 
0  1 


(5-29) 


In  view  of  the  definition  given  in  Equation  (5-5): 

nT 


XB  yB  zb]  = 


'"x  ",T 

XB 

yB 


-B 


=  rB] 


From  Equations  (5-29)  to  (5-31): 


where 


COS0Z 

-sin0z 

0‘ 

Rli 

R12 

R13 

sin0z 

cos0z 

0 

= 

R2i 

R22 

R23 

0 

0 

1 

R3i 

R32 

R33 

(5-30) 


(5-31) 


(5-32) 


(5-33) 
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Comparing  R®  in  Equation  (5-33)  with  C®  in  Equation  (5-6),  we  note  that  R®  is  equal  to  the 
transpose  of  C®: 


*i  = 

[<$f 

(5-34) 

c?  = 

Mf 

(5-35) 

From  Equation  (5-33): 

cos0z  -sin0z  0 
sin  0Z  cos  0Z  0 

0  0  1 

=  cos20z  +  sin20z  =  1  (5-36) 

It  follows  from  Equation  (5-8)  and  Equations  (5-34)  through  (5-36), 

|C|  =  |  CT|  =  |R|  =  |Rt|  =  1  (5-37) 

Perhaps  for  historical  or  implementation-related  reasons,  the  eye  movement  community 
seems  to  prefer  the  R  matrix  (described  in  Equations  (5-29)  through  (5-33))  corresponding  to  row 
vector  representation  of  coordinates,  while  virtually  all  other  scientific  and  engineering  communi¬ 
ties  prefer  the  C  matrix  (described  in  Equations  (5-4)  through  (5-7))  corresponding  to  column  vector 
representation  of  coordinates. 

5.3  BASIC  ROTATION  AROUND  THE  Y-AXIS 

Next,  we  consider  the  counterclockwise  rotation  around  the  common  Y-axis  of  Frame  B  rela¬ 
tive  to  Frame  A  by  0y. 

Following  exactly  the  same  procedure  as  before,  we  get 
xB  =  xAcos0Y  +  y  A0  —  zAsin0y 
yB  =  xAo  +  yA  +  zA0 

zB  =  xA  sin  0Y  +  yA0  +  zAcos0Y  (5-38) 
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or 


XB 

yB 

— 

zB 

cos0Y 
0 

sin  0v 


0 

1 

0 


—  sin  0y 
0 

cos0Y 


[ 


A 

Ya 

Zt 


It  follows,  using  notations  given  in  Equation  (5-5), 


rB  =  CBrA 


in  which  C?  is  given  by: 


COS0Y 

0 

-sin0Y" 

= 

0 

1 

0 

sin  0  Y 

0 

cos0Y 

CB 

'-'A 


for  Y-axis  rotation  by  0y. 

If  we  use  row  vector  representation,  we  get  from  Equation  (5-38): 

[xB  yB  zB]  =  [xA  yA  zA] 

or,  in  terms  of  notations  defined  in  Equation  (5-5): 

-iT  r  . tT 


COS0y 

0 

sin  0, 

0 

1 

0 

-  sin  0y 

0 

cos  0 

[rB]  =  [rA]  R 


where 


ra  = 


COS  0y 

0 

sin0. 

0 

1 

0 

-  sin  0y 

0 

COS0 

we  note,  comparing  Equation  (5-43)  with  Equation  (5-40): 


rb  =  [cb]  =  c 


(5-39) 


(5-40) 


(5-41) 


(5-42) 


(5-43) 


(5-44) 


as  before. 
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5.4  BASIC  ROTATION  AROUND  THE  X-AXIS 


Finally,  we  consider  the  counterclockwise  rotation  around  the  (originally)  common  X-axis 
of  Frame  B  relative  to  Frame  A  by  0X. 

Following  exactly  the  same  procedure  as  before,  we  get 

xb  =  xa  +  Ya  0  +  ZA° 

Yb  =  XA°  +  yAcos0x  +  zAsin0X 

zB  =  xA0  -  yAsin0x  +  zAcos0x  (5-45) 


From  Equation  (5-45), 


■kb' 

'l  0  O' 

■xa' 

Yb 

= 

0  cos  0X  sin  0X 

yA 

zb 

0  — sin  0X  cos  0X 

ZA 

(5-46) 


rB  =  cArA 


(5-47) 


We  may  also  express  Equation  (5-45),  using  row  vector  representation: 


[XB 


yb  zb]  =  [xa  yA 


0  0 
cos0x  -  sin  0X 
sin  0X  cos  0X 


(5-48) 


Since  [AB]T  is  equal  to  BTAT  (see  Appendix  A),  taking  the  transpose  of  Equation  (5-47): 


(5-49) 


Comparing  Equation  (5-49)  with  Equation  (5-48),  we  conclude 

[cyT  =  R5 


(5-50) 


as  shown  in  Equation  (5-34). 
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Taking  the  transpose  of  Equation  (5-50): 


T  T 

Since  [c®]  =  Cg  and  Cg  ]  =  C® ,  we  have: 


T 


(5-51) 


as  shown  in  Equation  (5-35). 

That  is,  C®  corresponding  to  the  column  vector  representation  is  the  transpose  of  R®  corre¬ 
sponding  to  the  row  vector  representation. 

5.5  COLUMN  VECTOR  REPRESENTATION  OF  TWO 
CONSECUTIVE  ROTATIONS 

Now  suppose  r®  in  the  Frame  B  is  transformed  to  rc  in  Frame  C  (not  to  be  confused  with 
Rotation  Matrix  C).  Then,  using  the  column  vector  representation: 


rc  =  Cg  r®  (5-52) 

Substituting  Equation  (5-47)  in  Equation  (5-52)  for  rB  : 
rc  =  Cg  [c»  rA] 

=  Cg  C\  rA  (5-53) 

It  is  also  true  that 

rc  =  Cg  rA  (5-54) 

From  Equations  (5-53)  and  (5-54),  we  get: 

Cg  =  CgC*  (5-55) 


Note  that  in  Equation  (5-55),  the  matrix  Cg  corresponding  to  the  second  rotation  is  multiplied 
from  the  left,  and  that,  while  the  rotations  are  sequential  and  might  appear  additive,  they  are  repre¬ 
sented  mathematically  by  multiplication  of  matrices.  Also,  note  that  sequence  in  which  the  rotations 
take  place  is  important  since,  in  general,  Cg  Cg  *  Cg  Cg. 
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5.6  ROW  VECTOR  REPRESENTATION  OF  TWO 
CONSECUTIVE  ROTATIONS 

Next,  using  the  row  vector  representation  for  the  transformation  from  Frame  B  to  Frame  C: 


[rc]T  =  [rB]T  R£  (5-56) 

Substituting  Equation  (5-48)  into  Equation  (5-56)  for  (rB)T : 

[rf  =  [[rf  R»]rC 

=  rA]T  Rb  Rg  (5-57) 

Since  it  is  also  true  that 
T  T 

[rc]  =  [rA]  Rj  (5-58) 

we  conclude  from  Equation  (5-57)  and  Equation  (5-58)  that  when  Frame  A  is  rotated  to  Frame  B, 
and  Frame  B  is  rotated  to  Frame  C,  the  rotations  from  Frame  A  to  Frame  C  equate  to: 

Ra  =  r!rb  (5-59) 


Note  that  in  Equation  (5-59),  the  matrix  Rg  for  the  second  rotation  based  on  the  transformation 

of  row  vectors  multiplies  from  the  right,  unlike  the  case  of  Cg  which  multiplies  from  the  left  as 
shown  in  Equation  (5-55). 

Using  similar  procedures,  we  can  show  that  when  Frame  C  is  rotated  to  Frame  D: 

C«  =  CgC^C»  (5-60) 

while 

RD  =  R?  rc  R»  (5-61) 
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Taking  the  transpose  of  Equation  (5-60),  and  using  the  formula  [AB]T  =  BTAT  (see 
Appendix  A), 

[c£f  =  [cgcfcjf 
=  [cJ]T  [eg  Cg]T 

=  [cjf  [cgf  [Cgf  (5-62) 

Using  Equation  (5-34)  and  its  extension  in  Equation  (5-62): 

R°  =  R®  R=  R°  (5-63) 

which  is  identical  with  Equation  (5-61). 
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SECTION  6 

ROTATION  MATRICES  FOR  FICK’S  SYSTEM  AND 
HELMHOLTZ’S  SYSTEM  OF  EYE  ROTATIONS 


In  the  following  pages,  we  will  derive  the  resultant  Rotation  Matrix  corresponding  to  the  three 
consecutive  Euler  angle  rotations  according  to  the  Fick’s  system  (discussed  in  Section  4.1).  We  will 
also  derive  the  resultant  Rotation  Matrix  corresponding  to  the  three  consecutive  Euler  angle  rota¬ 
tions  according  to  the  Helmholtz’s  system  (discussed  in  Section  4.2). 

We  assume  the  Head  Frame  and  the  Eye  Frame  are  initially  aligned  with  their  X-axes  pointed 
straight  forward,  their  Y-axes  horizontally  left,  and  their  Z-axes  vertically  upward. 

The  Rotation  Matrices  for  the  three  basic  rotations  around  the  Z-axis  (yaw  or  horizontal  rota¬ 
tion),  the  Y-axis  (pitch  or  vertical  rotation),  and  the  X-axis  (roll  or  torsional  rotation)  was  derived 
in  Section  5. 


6.1  ROTATION  MATRIX  FOR  FICK’S  SYSTEM 


In  Fick’s  system  of  Euler  Angles,  the  first  rotation  of  the  Eye  Frame  is  the  Zn-axis  of  the 
Head  Frame;  this  results  in  the  Frame  FI  (F  for  Fick).  The  second  rotation  of  the  Eye  Frame  is  about 
the  Ypi-axis  of  the  Frame  FI ;  this  results  in  Frame  F2.  The  third  rotation  of  the  Eye  Frame  is  about 
the  Xp2-axis  of  Frame  F2;  this  results  in  the  final  orientation  of  the  Eye  Frame  in  Fick’s  system. 

FICK  ROTATIONS:  (see  Section  5) 


1)  First  rotation  about  Zn-axis  by  67: 

Head  Frame  H  becomes  Frame  Fi  (F  for  Fick) 

rFl  _  rFl  rH 
r  ~  ''•'H  r 


Using  the  X-axis  rotation  matrix: 


rxpii 

COS0Z 

sin  0Z 

0 

hn  = 

-sin0z 

cos  0Z 

0 

M 

0 

0 

1 

2)  Second  rotation  about  Ypi-axis  by  0y: 
Frame  Fi  becomes  Frame  F2 

1 •F2  =  Cg  rF1 


(6-1) 


(6-2) 


(6-3) 
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Using  the  Y-axis  Rotation  Matrix: 

xF21  [cos0y  0  -  sin  0Y"|  rXpi- 

Yf2  =  0  1  0  Yfi  (6-4) 

Zp2_  sin0Y  0  cos0Y  _Zpl_ 

Substituting  Equation  (6-1)  into  Equation  (6-3): 

rF2  =  Cg  Cg  rH  (6-5) 

Notice  that  the  matrix  for  the  second  rotation  Cg  is  multiplying  the  matrix  for  the  first  rotation  Cg 
from  the  left. 

Substituting  from  Equation  (6-2)  and  Equation  (6-4)  into  Equation  (6-5): 

Xf2"|  [cos0Y  o  -sin0Yl  T  cos0z  sin0z  ol  rxH‘ 

yF2  =0  1  0  -sin0z  cos0z  0  yH  (6-6) 

Zp2_  sin0Y  0  cos0Y  0  0  1  .  H_ 

3)  Third  rotation  about  XF2-axis  by  0x: 

Frame  F2  becomes  Frame  F3  or  the  Eye’s  final  Frame  E. 
fE  _  j.F3  _  qF3  rF2 


Using  the  Z-axis  rotation  matrix: 


Substituting  Equation  (6-5)  into  Equation  (6-7)  for  rF2 : 

rE  =  eg  eg  eg  rH  (6-9) 


Notice  again  that  Cg  multiplies  CF2  from  the  left. 


6-2 
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Taking  the  transpose  of  Equation  (6-13) 


rEf  =  [cE  r")1 


•H]T[cgf 


rE]  =  [rH]  Rg 


(6-14) 


Equation  (6-14)  may  be  written  as 


XE  y£  ze]  =  [XH  yH  zh]  Rh 


(6-15) 


where  =  (C|)  is  the  transpose  of  the  Coefficient  Matrix  C|  of  Equation  (6-12)  for  the  Fick’s 
System,  which  is  given  below  in  Equation  (6-16). 


pE  _  FpiE 

KH  ~~  [UH 


(6-16) 


cos  0y  cos  0Z  -  cos  sin  +  sin  sin  ® y  cos  sin0x  sin  0Z  +  cos  0X  sin  0 Y  cos  0Z  ’ 

cos0Ysin0z  cos0xcos0z  4-  sin  0X  sin  0Y  sin  0Z  -sin0xcos0z  +  cos  0X  sin  0 Y  sin  0Z 
-sin0Y  sin0xcos0Y  cos0xcos0Y 


for  the  Fick’s  system  of  Euler  Angle  Rotations.  (This  is  the  transpose  (rows  and  columns  exchange) 
of  the  Coefficient  Matrix  of  Equation  (6-12).) 

Now,  taking  the  transpose  of  Equation  (6-9)  (see  Appendix  A): 

[rE]T  =  [cgc£CE'rH]T 

-  bfflGgC^CS1]7 

=  [«*1T  [cff  f  [eg  cgf 

=  [rH]T  [Ch  ]T  [cR]T  [cp|]T  (6-17) 
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or 


[rE]T  =  [rH]T[RpR£Rg] 


(6-18) 


or 


>F3 

"H 


[*T-[*T»1 

=  [>-h)T  Rf 

which  is  equal  to  Equation  (6-14)  and  Equation  (6-15). 


>E 

'‘H 


(6-19) 


6.2  ROTATION  MATRIX  FOR  HELMHOLTZ’S  SYSTEM 

In  Helmholtz’s  system  of  Euler  Angles,  the  first  rotation  of  the  Eye  Frame  is  around  the 
Yjj-axis  of  the  Head  Frame.  This  results  in  Frame  Hi  (H  for  Helmholtz).  The  second  rotation  of 
the  Eye  Frame  is  around  the  ZHi-axis  of  the  Frame  Hi.  This  results  in  the  Frame  H2.  The  third 
rotation  of  the  Eye  Frame  is  around  Xn2~axis  of  the  Frame  H2.  This  results  in  the  final  orientation 
of  the  Eye  Frame  in  Helmholtz’s  system. 

HELMHOLTZ  ROTATIONS:  (see  Section  5) 

1)  First  rotation  around  the  Yy-axis  by  0y : 

Frame  H  becomes  Frame  Hi  : 


rH1  =  Cg1  r 


,H 


(6-20) 


Using  the  Y-axis  rotation  matrix: 


XH1 

COS0y 

0 

Yhi 

— 

0 

1 

ZH1 

sin0Y 

0 

-sin0Y 

0 

COS0y 


2)  Second  rotation  around  Zm-axis  by  Qz'. 
Frame  Hi  becomes  Frame  H2 

_H2  _  rH2  -HI 
r  ~  '-'HI  r 


XH 

ZH 


(6-21) 


(6-22) 
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Using  the  Z-axis  rotation  matrix: 


XH2 

cos  0Z  sin  0Z  0 

XH1 

ym 

= 

-  sin  0Z  cos  07  0 

yHi 

ZH2 

3 

O 

1 

O 

ZH1 

Substituting  Equation  (6-20)  into  Equation  (6-22): 


„H2  _  CH2 


HI 


->H1  _H 
'H  r 


Substituting  Equation  (6-21)  and  Equation  (6-23)  into  Equation  (6-24): 


“ 

XH2 

ym 

= 

ZH2 

U  J 

_ 

cos  0Z  sin  0Z 
-  sin  0Z  cos  0Z 
0  0 


COS  0y 
0 

sin  0Y 


0 

1 

0 


-sin0Y 

0 

COS0V 


ah 

yH 

ZH 


3)  Third  rotation  around  Xn2-axis  by  0x; 

Frame  H2  becomes  Frame  H3  or  the  Eye’s  final  Frame  E. 

_E  _  _H3  _  rH3  _H2 
r  “  r  ~  UH2  r 

Using  the  X-axis  rotation  matrix: 


XE 

Ye 

= 

ZE 

■ 

0 

COS0-, 


0 

sin0x 
cos  0X 


XH2 

YHI 

ZH2 

Substituting  Equation  (6-24)  into  Equation  (6-26): 

_E  _  rH3  rH2  rHl  „H 
i  ^H2  ^ui  ^4 


"HI 


(6-23) 


(6-24) 


(6-25) 


(6-26) 


(6-27) 


(6-28) 
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SECTION  7 

ROTATION  MATRIX  FOR  SMALL  ANGLES 


Although  the  Rotation  Matrix  does  not  commute  for  finite  angles,  it  does  commute  for  very 
small  angles  for  which  the  approximations  sin  0  =  0,  cos  0  =  1,  and  00  — 0  are  justified.  Consider 

the  Rotation  Matrix  C|  (not  R^,  which  is  (c]|)  )  of  Fick’s  system.  Using  the  above  approxima¬ 
tions,  we  have  for  small  angle  rotations  of  Fick’s  system, 


-0Z  +  0X0Y 
®Y 


®z 

1  +  0X0Y0z 
-0X  +  0Y0Z 


(7-1) 


Since  for  very  small  angles  0x0y— 0  »  0x9y0z— 0  >  0x®z— 0  »  and  ©Y^Z— 0,  it  follows  from 
Equation  (7-1): 


1 

~®z 

0Y 


9Z  -0y 

1  0X 
-0X  1 


(7-2) 


for  Fick’s  system. 


T 

Now,  consider  the  Rotation  Matrix  (not  R^ ,  which  is  (C^)  )  of  Helmholtz’s  system  using 
the  approximations  sin  0  =  0,  cos  0  =  1,  and  00  =  0  for  small  angle  rotations: 


— 0Z  +  0X0Y 
9X®z  +  ®y 


1 

-0X 


-0Y 

0Z0Y  9 


X 


-0x0z0y  +  lj 


(7-3) 


for  Helmholtz’s  system.  Note  that  Equation  (7-3)  is  different  from  Equation  (7-1).  However,  with 
approximations  0x0y  =  0X0Z  =  0z0y  =  0x0y0z  =  0,  it  may  be  reduced  to: 


(7-4) 


which  is  identical  to  Equation  (7-2).  Equation  (7-4)  is  the  Rotation  Matrix  for  the  Helmholtz’s  sys¬ 
tem  with  angles  small  enough  to  justify  sin  0  =  0  and  cos0  =  1. 
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Thus,  for  small  angles, 

Cj|  (for  Fick’s  System)  =  (for  Helmholtz’s  System) 

1  0Z  -0Y 

=  ~ez  1  0X  (7-5) 

0Y  -0X  l 

That  is,  for  small  angles,  the  sequence  of  rotation  does  not  matter,  and;  therefore,  the  Rotation 
Matrices  commute,  or 

cg(9x)cR(eY)cS>(ez)  =  c»(ex)c“(ez)cH>(eY) 

making  the  final  eye  positions  equivalent  in  both  Helmholtz’s  and  Fick’s  systems  for  the  same 
amount  of  angular  rotations.  That  is,  for  the  separate  rotations  CVC2  and  C3  : 

CiC2C3  =  ci^3C2  =  C2CjC3  =  C2C3Cj  =  C3CjC2  =  C3C2Cj  (7-6) 

for  small  angle  rotations. 
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SECTION  8 

A  DIFFERENT  VIEW  OF  ROTATION  MATRIX 


Consider  a  vector  r  has  components  xa,  yA  and  za  in  Frame  A,  and  components  xb,  ys,  and 
zb  in  Frame  B.  Frame  B  is  displaced  counterclockwise  from  Frame  A  (reference  frame)  by  a  rota¬ 
tion  angle  0z  about  the  common  Z-axis.  For  this  case,  we  found  out  the  Rotation  Matrix  relating 
rA  =  [xa  yA  za]t  to  rB  =  [xb  yB  zb]t  is  given  by  (refer  to  Section  5): 


XB 

cos  0Z  sin  0Z  O' 

V 

yb 

= 

-  sin  0Z  cos  0Z  0 

yA 

ZB 

0  0  1 

ZA 

(8-1) 


Now  instead  of  rotating  Frame  B  relative  to  Frame  A,  a  vector  ri  in  Frame  A  with  coordinates 
xi,  yi,  and  zi  is  rotated  clockwise  about  the  negative  z-axis,  according  to  the  right  hand  rule  to  a 
new  position  r2  in  the  same  frame  with  coordinates  X2,  y2,  and  Z2  by  the  angle  0z-  This  means 
ri  is  moved  downward  on  the  X-Y  plane  as  shown  in  Figure  8-1. 


Referring  to  Figure  8-1,  0  is  the  angle  between  ri  and  T2-  (j>  is  the  angle  between  ri  and  the 
X-axis.  Note  that  ri  and  r2  have  the  same  magnitude,  which  we  denote  by  r. 

It  follows: 


xi  =  r  cos  <j) 

(8-2) 

yi  =  r  sin  <|> 

(8-3) 
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Using  the  trigonometry  formula  for  cos  (<}>— 0)  and  sin  (<f>-0): 

X2  =  r  cos  (<J>— 0) 

=  r(cos<))  cosG  +  sin<J)  sin0) 

=  cos0  (rcos<(>)  +  sin0  (rsin<}>) 

=  (cos0)xi  +  (sin0)yi 


Using  Equation  (8-2)  and  Equation  (8-3): 
y2  =  r  sin  (<})— 0) 

=  r(sin<(>  cos0  -  cos<J)  sin0) 

=  cos0  (rsin<j>)  -  sin0  (rcos<t>) 
=  (cos0)yi  -  (sin0)xi 
=  -(sin0)xi  +  (cos0)yi 


Zi  =  z2 
It  follows: 


"x2* 

cos  0  sin  0  0* 

V 

yi 

= 

-sin  0  cos0  0 

yi 

z2 

0  0  1 

zi 

(8-4) 


(8-5) 

(8-6) 


(8-7) 


Comparing  the  Coefficient  Matrix  of  Equation  (8-7)  with  that  of  Equation  (8-1),  we  notice  that 
they  are  identical. 

Generalizing  the  above  results,  the  Rotation  Matrix  corresponding  to  the  rotation  of  one  frame 
relative  to  the  reference  frame  is  the  same  as  the  Rotation  Matrix  corresponding  to  the  rotation  of 
a  vector  in  the  opposite  direction  (with  same  angle  of  rotation)  within  the  original  reference  frame. 
Most  of  the  time,  the  Rotation  Matrix  is  used  in  the  former  context.  However,  if  it  is  used  in  the  latter 
context,  it  should  be  made  explicit  to  avoid  confusion. 
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SECTION  9 

DETERMINATION  OF  ROTATION  ANGLES  FROM 
ROTATION  MATRIX  OF  FICK’S  SYSTEM 


Until  now  we  have  defined  angular  rotations  and,  from  these  rotations,  have  determined  the 
elements  of  a  Rotation  Matrix.  Now  we  reverse  the  process,  and  assume  we  have  a  Rotation  Matrix 
and  wish  to  obtain  the  rotation  angles  that  correspond  to  the  matrix. 


9.1  ANGLES  OF  FICK’S  SYSTEM  MATRIX 

(SOME  OF  TfflS  SECTION  HAS  BEEN  DISCUSSED  IN  AN  EARLIER  SECTION) 


Assume  the  Eye  Frame  is  rotated  to  another  orientation  by  rotating  the  initial  orientation 
around  its  X-axis  by  angle  T  (for  Torsion).  The  corresponding  Rotation  Matrix  CX(T)  is  given 

by: 


CX(T)  = 


1  0  0 
0  cos  T  Sin  T 
0  -sin  T  cos  T 


(9-1) 


which  is  Equation  (5-46). 

The  Rotation  Matrix  Cy  (V)  around  the  Y-axis  by  angle  V  (for  vertical  motion  of  the  Eye-axis) 
is  given  by: 


Cy  (V) 


cos  V  0  -sin  V 

0  1  0 

sinV  o  cosV 


(9-2) 


which  is  Equation  (5-39). 


The  Rotation  Matrix  Cz  (H)  around  the  Z-axis  by  angle  H  (for  horizontal  motion  of  the  Eye- 
axis)  is  given  by: 


CZ(H) 


cos  H  sin  H  0 

-sin  H  cos  H  0 

0  0  1 


(9-3) 


which  is  Equation  (5-4). 
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The  Equation  (9-7)  is  the  final  Rotation  Matrix  of  Fick’s  system.  Equation  (9-7)  is  the  same  as  the 
Coefficient  Matrix  derived  in  Equation  (6-12).  From  Equation  (9-7): 

C13  =  —  sinV  =  sin(—  V) 


or 


sin  1  C13  =  -  V 


or 


V  =  -  sin-1  C13 


(9-8) 


Thus,  we  can  determine  the  value  of  V  from  the  value  of  C13  determined  from  experiments. 
Next,  from  Equation  (9-7): 

C12  =  cosV  sin  FI 


sinFl  = 


C12 
cos  V 


or 


H 


(9-9) 


where  cos  V  is  determined  using  V  given  in  Equation  (9-8).  Note  that,  both  in  Fick’s  system  and 
Helmholtz’s  system,  all  angles  of  eye  rotations  are  physiologically  constrained  to  be  much  less  than 
90°.  Therefore,  cos  V  will  never  be  0. 

Now,  from  Equation  (9-7): 

C33  =  cosT  cosV 


or 


or 


cosT  = 


C33 

cosV 


T 


=  cos 


C33 

cosV 


(9-10) 
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The  Torsion  angle  T  may  also  be  determined  from  Equation  (9-7): 


C23  =  sinT  cosV 


sinT  = 


'23 


cos  V 


qp  -  - 1  ^23 

T  =  sin  1  - ^7 

cos  V 


Now  we  know  C  =  RT . 
That  is: 


C11 

C12 

C13] 

[*11 

R12 

Rnl 

C-21 

^22 

^23 

R21 

R22 

R23 

^31 

^32 

C33^ 

R31 

R32 

R33 

From  (9-12),  we  get: 

Cj3  =  R31;  C12  =  R2i;  C23  =  R32;  C33  =  R33 
From  Equation  (9-9)  and  Equation  (9-13), 

H  =  sin-'(^|V) 

Similarly, 

T  —  cos  (  cosv  )  or  T  —  sin  (  cos  v  ) 

V  =  -sin-1  R31 


(9-11) 


(9-12) 

(9-13) 

(9-14) 

(9-15) 

(9-16) 
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SECTION  10 

DETERMINATION  OF  ROTATION  ANGLES  FROM 
ROTATION  MATRIX  OF  HELMHOLTZ’S  SYSTEM 


10.1  ANGLES  OF  HELMHOLTZ’S  SYSTEM  MATRIX 


By  definition,  in  the  Helmholtz’s  System,  the  sequence  of  rotation  is:  Y-axis  rotation,  fol¬ 
lowed  by  Z-axis  rotation,  followed  by  X-axis  rotation.  Thus,  the  final  Rotation  Matrix  C  in  the 
Helmholtz’s  System  is  given  by,  using  the  notations  described  in  the  previous  section: 


C  =  CX(T)CZ(H)CY(V) 


Now:  . 


cosH 

sinH  0 

cosV  0  -sin  V 

Cz(H)Cy(V)  = 

-sin  H 

0 

cosH  0 

0  1 

0  1  0 

sin  V0  cos  V 

cos  H  cos  V 
-  sin  H  cos  V 
sin  V 


sinH  -cosHsinV 
cosH  sin  H  sin  V 

0  cos  V 


(10-1) 


(10-2) 


(10-3) 


10-1 


It  follows: 

1  0  o  1  cos  H  cos  V  sin  H  -cos  H  sin  V 

Cx(t)Cz(h)Cy(v)  =  0  cosT  sinT  -sin  H  cos  Y  cosH  sin  H  sin  V  (10-4) 

0  -sinT  cosT  sinV  0  cosV 
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where  cos  H  is  determined  using  H  in  Equation  (10-6).  Note  that,  both  in  Fick’s  System  and  Helmholtz’s  System,  all  angles 
of  eye  rotations  are  constrained  to  be  less  than  90°,  so  that  cos  H  will  never  be  0. 
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Next, 


C22  =  cos  T  cos  H 


cosT  = 


_  *-22 


cosH 


T  =  cos  1 


-22 


cosH 


Torsion  T  may  also  be  found  from: 


-32 


=  -sinT  cos  H 


C32  =  sin  (-T)  cos  H 


sin  (-T)  = 


-32 


cosH 


<-T)  ■ si"-'  tdk 


=  -sin 


_1  C32 


cosH 


Now,  we  know: 


fRll 

R12 

R13 

C11 

C21 

C3l] 

R21 

R22 

R23 

*-12 

*-22 

*-32 

R31 

R32 

R33 

C13 

C23 

C33 

From  (10-10): 

Cl2  =  R21’  ^13  =  R31’  ^22  =  R22’  ^32  =  R23 

It  follows  from  (10-6)  to  (10-9),  using  (10-11): 

H  =  sin-1  Ro 


V  =  -sin 


'■21 


-1  R31 


cosH 


T  =  cos  1  — or  T  =  -sin  1  —^77  ■ 
cosH  cosH 


(10-8) 


(10-9) 

(10-10) 

(10-11) 

(10-12) 

(10-13) 

(10-14) 
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SECTION  11 

ORTHOGONAL  MATRIX  AND  ORTHOGONALITY  CONDITION 


In  Section  5,  we  have  shown  by  means  of  a  simple  example  based  on  one  axis  rotation  that  the 
Coordinate  Transformation  Matrix  (Rotation  Matrix)  from  Frame  A  to  Frame  B,  denoted  by 
C® ,  has  a  property  that  its  inverse  is  equal  to  its  transpose.  That  is, 

(C =  (cl)\ 

where  C®  is  an  Orthogonal  Matrix. 

It  follows: 

(C  w1  =  C  J  (c*)T 

Also, 

nr'ci  =  (c»)Tcs  = 

In  this  Section,  we  will  show  that  Equation  (11-1)  is  valid  for  any  arbitrary  orientation  of  one 
frame  relative  to  another  frame,  regardless  of  how  frames  may  have  reached  the  mutual  orientation. 

Consider  two  frames.  Frame  A  and  Frame  B.  Frame  A  has  Xa,  Ya,  and  Za  axes  with  unit 
vectors  I,  J,  and  K,  respectively.  Frame  B  has  Xg,  Yg,  and  Zg  axes  with  unit  vectors  i,  j,  and  k, 
respectively.  An  arbitrary  vector  has  components  xa,  yA>  and  ZAin  Frame  A,  and  components  xg, 
yg,  and  zg  in  Frame  B  as  shown  in  Figure  11-1. 


_  pB  p^A  _  p^B  _  t 


A  p^B 
^B  UA 


cj  =  i 


(11-2) 


(11-3) 


11-1 
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Figure  11-1.  An  Arbitrary  Vector  with  Components  Expressed  in  Frames  A  and  B 


We  want  to  determine  xB,  yB  and  zB  of  Frame  B  in  terms  of  xA,  yA,  and  zA  of  Frame  A.  First 
consider  xB,  which  consists  of, 

xB  =  (component  of  xA  along  i  or  XB  axis) 

+  (component  of  yA  along  i  or  XB  axis) 

+  (component  of  zA  along  i  or  XB  axis) 
xB  =  xAcos  (i,I)  +  yA  cos(i,J)  +  zAcos  (i,K)  (11-4) 

in  which  cos  (i,  I)  is  the  cosine  of  angle  formed  by  unit  vector  i  and  I,  and  so  on. 

Similarly, 

yB  =  xAcos  (j,I)  +  yAcos  (j,J)  +  zAcos  (j,K)  (11-5) 

zB  =  xAcos  (k,I)  +  yAcos  (k,J)  +  zAcos  (k,K)  (11-6) 


11-2 
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From  Equation  (11-4)  to  Equation  (11-6) 

cos(i,I)  cos(i,J)  cos(i,  K)'B  xa' 

cos(j,  I)  cos(j,  J)  cos(j,K)  yA  (H-7) 

cos(k,I)  cos(k,J)  cos(k,K)  z* 

J  l  J 

The  Coefficient  Matrix  of  Equation  (11-7)  above  is  called  Direction  Cosine  Matrix  (which 
is  also  called  Rotation  Matrix  or  Coordinate  Transformation  Matrix)  because  all  of  its  elements 
are  cosines  of  the  angle  between  directions  of  axes  of  two  different  frames. 

By  definition  as  we  discussed  in  Section  1, 

i  •  I  =  |i|  |I|  cos  (i,I)  =  cos  (i,  I) 
i  •  J  =  |i|  |J|  cos  (i,J)  =  cos  (i,J),etc. 
since  i,  j,  k  and  I,  J,  K  are  unit  vectors. 

It  follows  that  Equation  (11-7)  may  be  written  as 

i  •  I  i  J  i  •  K]B  fxA 

J  *  1  J  *  J  J  •  K  Ya  (11-8) 

k  I  k  •  J  k  K  zA 

Ja  l 


(11-9) 


in  which  the  meaning  of  Cn,  C12,  etc.  are  obvious  by  comparing  Equation  (11-7)  and  Equa¬ 
tion  (11-8)  to  Equation  (11-9). 

We  may  summarize  Equation  (11-8)  and  Equation  (11-9)  by 

xb' 

yB  =C« 
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or 

rB  =  CB  rA 

where  CB  represents  the  Coefficient  Matrix  on  the  right  side  of  Equation  (11-9). 
Next,  we  reverse  the  process.  We  want  to  find  xA,  yA,  zA  in  terms  of  xB,  yB,  zB. 
xA  =  (component  of  xB  along  I  or  XA  axis) 

+  (component  of  yB  along  I  or  XA  axis) 

+  (component  of  zB  along  I  or  XA  axis) 

xA  =  cos  (i,  I)  xB  +  cos  (j,I)  yB  +  cos  (k,I)  zB 


Similarly, 

yA  =  cos  (i,  J)  xB  +  cos  (j,  J)  yB  +  cos  (k,  J)  zB 
zA  =  cos  (i,K)  xB  4-  cos  (j,K)  yB  +  COS  (k,K)  zB 


From  Equation  (11-12),  Equation  (11-13),  and  Equation  (11-14): 


V 

r  cos  (i,  I) 

cos(j,  I) 

cos  (k,  I) ' 

A 

V 

Ya 

cos(i,  J) 

cos(j,  J) 

cos  (k,  J) 

Yb 

ZA 

cos(i,K) 

cos(J,  K) 

cos  (k,  K) 

> 

B 

ZB 

i  I 

j  I 

k  r 

A 

V 

i  •  J 

j  J 

k  J 

Yb 

i  K 

j  K 

k  •  K 

D 

ZB 

It  follows 


V 

XB 

Ya 

<P3 

U 

II 

Yb 

ZA 

ZB 

or 


rA  =  CB  rB 


(11-10) 


(11-11) 

(11-12) 

(11-13) 

(11-14) 

(11-15) 

(11-16) 


(11-17) 
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where  Cg  represents  the  Coefficient  Matrix  on  the  right  sides  of  Equation  (11-15)  and  Equa¬ 
tion  (11-16). 


Comparing  Cg  given  in  (11-15),  with  C®  given  in  (11-8),  we  note  that  Cg  is  equal  to  the 
Transpose  (rows  and  columns  exchanged)  of  C®  .  That  is, 


Cg  =  (C?)T  and  (C®)  =  (Cg)T. 

(11-18) 

Substituting  Equation  (11-10)  into  Equation  (11-17): 

rA  =  c£  r® 

=  Cg  C®  rA. 

(11-19) 

Since  rA  =  I  rA  where  I  =  Identity  Matrix,  we  conclude: 

pA  pB  _  t 
'“'B  '-'A  l- 

(11-20) 

Using  Equation  (11-18)  in  Equation  (11-20): 

(C®)T  (C®)  =  I. 

(11-21) 

Since  (C®)"1  (C®)  =  I, 

(11-22) 

we  conclude: 

(C®)-1  =  (C®)T  . 

(11-23) 

Since  Equation  (11-23)  is  true  for  any  two  Orthogonal  Rotation  Matrixes,  we  have,  dropping 
subscripts  and  superscripts:  C-1  =  CT . 

It  follows: 

C^C  =  CTC  =  I 

(11-24) 

0 

0 

II 

o 

0 

H 

II 

hH 

(11-25) 

Although  Equation  (11-18)  through  Equation  (11-23)  were  shown  in  Section  5  using  a  single 
axis  rotation,  we  are  confirming  them  here  for  general  case. 

Cji  Cj2  C13 

Thus,  for  C  =  C21  C22  C23  (11-26) 

C31  C32  C33 

k.  J 
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we  have,  using  Equation  (11-24): 


C11 

C-21 

C31 

^12 

C13] 

’  1 

0 

°1 

CTC  = 

C12 

C22 

C32 

C21 

C22 

C23 

- 

0 

1 

0 

Cj3 

^23 

C33 

C31 

^32 

C33 

> 

0 

0 

1 

(11-27) 


Carrying  out  matrix  multiplication: 


CTC  = 


C?1  +  C21  +  C31 

Ci2C„  +  C22C21  +  c32c31 

^13^11  +  ^23^21  +  ^33^-31 


cnc12  +  c 


c„  +  c„c 


21^22 


31^32 


C11C13  +  C21C23 


+  c31c33 


n2  4.  f2  4.  p 2 

'-'12  ~  *~22  T  '-'32 
^13^12  +  ^23^22  +  ^-'33^32 


^12^13  +  ^22^23  C32C33 

^13  ■*"  ^23  *-% 


r  1  0 

0  1 

0  0 


(11-28) 


Equating  (1,1)  elements: 

Cfi  +  C221  +  C231  =  1  (11-29) 

Equating  (2,2)  elements: 

C12  +  C22  +  C32  =  1  (H-30) 

Equating  (3,3)  elements: 

CI3  +  C23  +  C^3  =  1  (11-31) 

Equating  (1,2)  elements: 

CnCj2  +  C21C22  +  C31C32  =  0  (11-32) 

Equating  (1,3)  elements: 

CijCis  +  C21C23  +  C31C33  =  0  (11-33) 

Equating  (2,3)  elements 

C12C13  +  C22C23  +  C32C33  =  0  (11-34) 


Referring  to  the  matrix  of  CTC  given  in  Equation  (11-28),  we  note  that  its  (2,1)  element  is  iden¬ 
tical  to  its  (1,2)  element;  (3,2)  element  to  (2,3)  element;  and  (3,1)  element  to  (1,3)  element. 
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If  we  use  Equation  (11-25),  we  have 


C11 

C-12 

( 

r>  " 

^13 

^11 

^21 

C31I 

CCT  = 

C21 

C22 

< 

^23 

C12 

C22 

C32 

C31 

k 

^32 

< 

^33^ 

^13 

C23 

C33 

'  1 

0 

01 

0 

1 

0 

0 

0 

1 

(11-35) 


which  leads  to  the  following  relationships  among  the  coefficients,  sometimes  called  the  redundan¬ 
cy  relationships. 

It  follows  from  Equation  (11-35): 


C2  +  c2  +  C2  =1 
Ml  T  M2  M3  1 


^21  +  ^22  +  ^23  ~  1 
C31  +  C32  +  C33  =  1 


(11-36) 

(11-37) 

(11-38) 


Equating  (1,2)  elements: 


C  11^21  ^-'12^22  C13C23  —  0 


(11-39) 


Equating  (2,  3)  elements: 


^21^31  ^22^32  ^23^33  —  ®  (11-40) 

Equating  (1,3)  elements: 

C„C31  +  C12C32  +  C13c33  =  0  (11-41) 

Equation  (11-29)  through  Equation  (11-34)  or,  equivalently.  Equation  (11-36)  through 
Equation  (11-41)  are  results  obtained  from  Equation  (11-28)  or  Equation  (11-35),  repeated  below: 

CTC  =  CCT  =  I  (11-42) 


which  is  called  orthogonality  condition. 

When  the  Rotation  Matrices  are  obtained  experimentally  or  determined  computationally  with 
computation  errors,  they  generally  do  not  satisfy  Equation  (11-42).  In  those  cases,  Equation  (11-29) 
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through  Equation  (11-34)  or  Equation  (11-36)  through  Equation  (11-41)  are  used  to  make  the 
Rotation  Matrices  satisfy  Equation  (11-42).  That  is,  the  equations  “normalize”  and  “orthogonalize” 
the  Rotation  Matrices. 

The  3x3  Direction  Cosine  Matrix  C  has  nine  direction  cosines.  Six  orthogonality  relations 
—  Equation  (11-29)  through  Equation  (11-34)  or  Equation  (11-38)  through  Equation  (11-41),  — 
contain  all  nine  elements.  Therefore,  these  six  equations  may  be  solved  for  six  elements  in  terms 
of  three  remaining  elements.  This  means  there  are  only  three  independent  parameters.  Since  we 
cannot  reduce  the  number  of  the  independent  parameters  to  less  than  three,  making  the  number  3 
the  minimum  number  of  independent  parameters  to  specify  the  rotation,  a  set  of  three  coordinate 
is  called  generalized  coordinates.  Although  there  are  a  number  of  such  sets  of  parameters  available, 
the  most  widely  accepted  are  Euler  Angles. 

A  more  useful  set  of  relationship  than  those  given  by  Equations  (11-39)  through  (11-41)  fol¬ 
lows  from  the  orthogonality  of  the  unit  triads  (i,  j,  k)  and  (I,  J,  K).  By  definition, 

i=jxk;  j  =  k  x  i;  k  =  ixj  (11-43) 


Applying  Equation  (11-9)  between  (i,  j,  k)  and  (I,  J,  K), 


f  .  5 

1 

Cjj  C12  Cj3 

r 

j 

= 

^21  C22  C-23 

J 

k 

C31  c32  c33 

K 

It  follows  from  Equation  (11-44): 

i  =  CnI  +  C12J  +  CjjK 
j  =  C2il  +  C22J  +  C23K 


k  =  C31I  +  C32J  +  C33K 


(11-44) 


(11-45) 


Using  Equation  (11-45)  in  Equation  (11-43)  for  i  =  j  X  k  : 

Cjj  I  +  C12  J  +  C13  K 

=  (c21 1  +  C22  J  +  C23  K)  x  (C31 1  +  C32  J  +  c33  k) 


C2 

c. 


J 

K 

^22 

C23 

-  (c22c33  C32  C23)  I 

^32 

c33 

— 

C21  ^33)  J  +  (C21  C33  C22  ^3l)  K 

(11-46) 

11-8 


NSW  CDD/MP-99/17 


Equating  coefficients  for  I,  J,  K  on  both  sides  of  Equation  (11-46),  we  get: 
Cn  =  C22C33  -  C32C23 
Cj.2  =  ^23^31  -  ^21^33 
C-13  =  ^-21^32  -  C3iC22- 

Similarly  we  obtain  from  j  =  k  x  i: 

C21  =  c32c13  -  C12C33 

^22  =  C11C33  -  C3iC13 
C23  =  C31C12  ~  C11^32- 

and  from  k  =  i  x  j  : 

C31  =  ^12^23  —  ^-22^13 
^32  =  ^21^13  -  ^11^23 
^33  =  ^11^22  -  ^21^12- 

Now  going  back  to  Equation  (11-32)  and  solving  for  Ch: 
n  _  ~  (^21^22  +  C31C32) 


(11-47) 


(11-48) 


(11-49) 


(11-50) 


If  C12  is  0,  Cn  cannot  be  obtained  from  Equation  (11-50).  However,  it  can  be  obtained  from 
the  first  equation  (Equation  (11-47))  without  these  restrictions.  In  this  sense,  Equations  (11-47), 
(11-48),  and  (11-49)  are  more  useful  than  Equations  (11-32),  (11-33),  and  (11-34)  because  they 
allow  us  to  avoid  singularities. 
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SECTION  12 

NORMALIZATION  OF  VECTOR  AND  ROTATION  MATRIX 


12.1  NORMALIZATION  OF  VECTORS 

Given  a  vector  r, 

r  =  ix  +  jy  +  kz 

we  want  to  normalize  r  so  that  its  magnitude  |r|  =  1 
That  is,  we  require 

M  =  I**  +  jy  +  M 

=  (x2  +  y2  +  z2)1/2  =  i 
Now  define  x,  y  ,  z  ,  and  r  by: 


r  =  ix  +  jy'  +  kz’ 


Then 


(x')2  +  (y)2  +  (z)2  =  4  +  —1  +  4 

|r|2  Irl2  Irl2 


x2  +  y2  +  z2 


_  x2  +  y2  +  z2 
x2  +  y2  +  z2 


=  1 


It  follows  that 


r  = 


(x'f  +  (y)2  +  (z'f  -  1 


which  implies  that  r  is  normalized. 


(12-1) 


(12-2) 


(12-3) 


(12-4) 
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12.2  NORMALIZATION  OF  MATRICES 
For  an  Orthogonal  Matrix: 

C11  Cl2  ^13 

c  =  C21  c22  c23 

C31  ^32  C33 

such  that 


C11 

C21 

C31 

R11 

R12 

Rj3 

CT  = 

C12 

C22 

C32 

— 

R21 

r22 

R23 

C13 

^23 

C33 

R31 

R32 

R33 

we  know  if  Equation  (12-5)  is  errorless: 

CCT  =  I  and  CTC  =  I 


(12-5) 


(12-6) 


Cn  C12  c13  Cn  C21  C31  |"l  0  0 

ccT  =  C21  c22  C23  C12  C22  C32  =  0  1  0=1 

C31  c32  c33  c13  c23  c33  0  0  1 

It  follows  from  the  above  equation: 

^11  +  ^12  +  ^13  =  * 

^21  ■*"  ^22  ^23  =  ^ 

C?1  +  C?2  +  c*3  =  1 


(12-7) 


(12-8) 


Equation  (12-8)  is  called  Normalization  condition. 

In  some  cases,  a  Rotation  Matrix  may  not  obey  Equation  (12-6).  The  next  paragraphs  show 
how  to  ensure  that  each  diagonal  element  of  CCT  is  equal  to  I. 

Now,  assume  the  Qj’s  in  Equation  (12-5)  have  errors. 

We  want  to  normalize  the  C  matrix  so  that  Equation  (12-8)  is  satisfied. 
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Now  define  Cri,  Cr2  and  Cr3  by: 

Cri  —  [^11  +  ^12  +  ^13] 
Cr2  —  [^21  +  ^"22  + 

Cr3  —  [^"31  +  C32  +  C33] 
and  also  define  Cn,Ci2  and  C13  by: 


C'  A 

^  11  =  c 


R1 


r>  A  Ci2 
12  =  c 

^R1 

r'  A  Ci3 

u  13  =  r 

'-'Rl 


where  the  symbol  A  denotes  the  definition. 


It  follows: 


,  ,  ,  c?,  +  c?~  +  c?_ 

[Cn]2  +  [C12]2  +  [C13]2  =  11  122 - ^  =  C21  +  C2  +  C23 

MRl  11  x  ^12  ^13 


Define  C2i,C22  and  C23  by: 


c2IA^ 

CR2 

ci2A£s 

^R2 


c23A^ 

'~R2 


(12-9) 


(12-10) 


1  (12-11) 


(12-12) 
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It  follows: 


[C2iV  +  [C22V  +  [C23V  = 


p2  +  r2  +  ri 
2  _  21  22  23 


[CR2]2 

Cll  +  +  C23 

C21  +  C22  +  C23 


=  1 


t  _ t 


Finally,  define  C31,  C32  and  C33  by: 


,  _  C31  ,  _  C32  ,  _  C33 

^31  p  ?  ^32  p  »  ^33  p 

^R3  CR3  ^R3 


It  follows: 


[c3i]2  +  [c32]2  +  [c33]2 


C31  +  C32  +  C33 

[CR3]2 

r2  4.  r2  4.  r2 
^31  32  33 

r2  4.  r2  4.  r2 
31  32  33 


1 


Therefore,  the  new  matrix  C  given  by: 


C'n 

C  12 

c'13' 

c'  = 

C21 

C'22 

c'23 

C31 

c'32 

c'33 

(12-13) 


(12-14) 


(12-15) 


(12-16) 


will  satisfy  Equation  (12-8)  with  Qj  replaced  by  C  y,  which  means  C  is  normalized. 
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SECTION  13 

ORTHOGONALIZATION  OF  VECTOR  AND  ROTATION  MATRIX 


13.1  ORTHOGONALIZATION  OF  VECTORS 

Given  a  set  of  three  non-orthogonal  vectors  with  a  common  origin,  we  want  to  find  a  set  of 
three  mutually  orthogonal  vectors  to  form  an  orthogonal  frame  of  reference. 

That  is,  from  three  non-orthogonal  axes  Xi,  X2  and  X3,  we  want  to  find  three  orthogonal  axes 
yi,  y2  and  y3.  Figure  13-1  shows  xi  and  X2,  which  are  not  perpendicular  to  each  other. 


P 


Figure  13-1.  Two  Non-Orthogonal  Vectors  xi  and  X2. 

Initially  we  choose  yi  to  be  xi,  or  xi  =  yi.  Then  we  draw  a  perpendicular  line  from  p,  (the  tip 
of  X2),  to  the  Xi  line,  and  call  the  intersection  of  this  line  with  xi  line  p .  Let 
Op'  =  yi,  and  p'p  =  y2 . 


Define  ki  by 

yi  =  kiyi 

where  kj  is  a  scaling  factor  for  yj. 
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It  follows  from  Figure  13-1,  using  vector  addition: 

X2  =  yi  +  y2  =  kiYi  +  y2 


(13-1) 


From  Equation  (13-1): 


y2  =  x2  -  kiyi 


(13-2) 


Since  yi  is  perpendicular  to  y2,  we  require  that  yi  •  y2  =  0 
or,  using  Equation  (13-2): 


yi  *  y2  =  yi  •  (*2  -  kiyi)  =  o 


(13-3) 


It  follows: 

yi  •  y2  =  yi  •  x2  -  kiyi  •  yi  =  0 
which  results  in: 

_  yi  •  x2 
kl  "  yrTT 


(13-4) 


(13-5) 


Substituting  Equation  (13-5)  into  Equation  (13-2)  for  kj: 

yi  •  x2 

y2  =  x2  -  y-^yi 


(13-6) 


•  As  a  check  to  see  if  yi  is  perpendicular  to  y2: 

(  yi  •  x2  \ 

yi-y2  =  yf(x2-yp77yij 

yi  •  x2 

=  yi  •  x2  -  y-Ty^yi  •  yi 

=  yi  •  x2  -  yi  •  x2  =  0 

which  shows  that  y2  as  determined  by  Equation  (13-6)  is  indeed  perpendicular  to  yi.  Now  that  we 
found  y2  determined  by  Equation  (13-6)  is  perpendicular  to  yi,  our  next  task  is  to  find  which  is 
perpendicular  to  both  yi  and  y2. 

Referring  to  Figure  13-2,  draw  a  perpendicular  line  from  Q  (the  tip  of  X3)  to  the  yj  -  y2  plane, 
and  call  the  intersection  Q  .  Drop  a  perpendicular  line  from  Q  to  Oy  j  line  and  call  the  intersection 
yi .  Finally,  draw  a  perpendicular  line  from  Q  to  Oy2  line  and  call  the  intersection  y2. 


(13-7) 


(13-8) 
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Figure  13-2.  Orthogonal  Set  of  Vectors 


We  define  k2  and  k3  by 

yi  =  k2yj  (13-9) 

and 

72  =  k3y2  (13-10) 

where  k2  and  k3  are  scaling  factors  for  yj  andy2  respectively. 

By  simple  vector  addition: 

OQ  =  yi  +  72  =  k2yi  +  k3y2.  (13-11) 

It  follows: 

x3  =  OQ'  +  QQ  =  k2y}  +  k3y2  +  y3.  (13-12) 

From  Equation  (13-12): 

y3  =  x3  -  k2yj  -  k3y2.  (13-13) 
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We  want  y2  and  y3  to  be  perpendicular  to  each  other,  or  y2  •  y3  =  0,  or  using  Equa¬ 
tion  (13-13): 


y2  •  y3  =  y 2  •  (x3  -  k2yi  -  k3y2) 

=  y2  •  x3  -  k2y2  •  yl  -  k3y2  •  y2  =  0 


(13-14) 


Since  we  have  shown  previously  that  y2  •  yi  =  0,  it  follows  from  Equation  (13-14): 

k3y2  •  y2  =  y2  •  x3  (13-15) 


or 


k,  = 


y2  *  x3 


(13-16) 


or 


3  y2  •  y2 

Finally,  we  want  yj  •  y3  =  0. 

Using  Equation  (13-13)  for  y3: 

yi  ■  y3  =  yi  •  (x3  _  k2yi  -  k3y2)  (13-17) 

=  yi  •  x3  -  k2yj  •  yj  -  k3yi  •  y2 

=  0  (13-18) 

Again,  since  we  shown  previously  that  yj  •  y2  =  0,  we  get  from  Equation  (13-18): 

yi  *  x3  —  k2yi  •  yi  =  0  (13-19) 

yi  •  x3 

k2  =  yi  •  yj  (13-20) 
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As  a  check  to  see  if  y3  is  perpendicular  to  yl5  using  Equation  (13-13)  for  y3  and  Equa¬ 
tion  (13-20)  for  k2: 

yi  •  y3  =  yi  •  (x3  -  k2yi  -  k3y2) 

=  yi  •  x3  -  k3yi  *  yi  -  k3yi  •  y2 


=  yrx3~  k2yi  *  yi 


yi  ‘x3 

=  yi‘x3~  yi  *  yx  yi  ’  yi 

=  yi  •  X3  -  yi  •  x3  =  0  (13-21) 

as  expected  because  yj  *  y2  =  0. 

To  see  if  y3  is  perpendicular  to  y2  using  Equation  (13-13)  for  y3  and  Equation  (13-16)  for  k3: 


y2  •  y3  =  y2  •  (x3  -  k2yi  -  k3y2) 


=  y2  •  x3  ~  hY2  •  yi  -  k3y2  *  y2 


=  y2  •  x3  -  k3y2  •  y2 


y2  *  x3 

=  y2  *  x3  -  y2  •  y2  y2  •  y2 

=  y2  •  x3  -  y2  •  x3  =  0  (13-22) 

as  expected  because  y2  *  yj  =  y^  •  y2  =  0. 

In  summary,  given  a  set  of  non-orthogonal  vectors  Xi,  X2,  x3,  we  can  find  the  corresponding 
orthogonal  set  yi,  y2,y3.  First,  we  let  yi  =  xi.  Then 


y2  =  x2  -  kiyi 


y3  =  x3  -  k2yr  -  k3y2 
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13.2  ORTHOGONALIZATION  OF  ROTATION  MATRICES 


Referring  to  Section  11,  consider  matrix  C: 

=  (Cj  c2c3) 


where 


;11 

Cl2 

Cj3 

21 

C22 

C23 

31 

C32 

C33 

column  1  :  Cj  = 


'll 

'21 

'31 


(13-23) 


column  2  :  Co  = 


'12 

'22 

'32 


column  3  :  C3  = 


'13 

-23 

'33 


It  follows: 


The  orthogonality  condition  states: 
CTC  =  I 


or 


C„ 

C21 

C31 

C11 

Cj2 

Cl3] 

'  1 

0 

0] 

C12 

C22 

C32 

C21 

C22 

C23 

=r 

0 

1 

0 

Cl3 

C23 

C33 

C31 

C32 

C33 

- 

0 

0 

1 

Using  Equation  (13-24)  and  Equation  (13-25)  in  Equation  (13-27): 


CTC  = 


CV 

Cl 

Cj 


(Cl  c2  c3) 


(13-24) 


Cl  -  (cn  C21  c31) ;  C2  -  (C12  C22  C32  ) ;  C3  -  (C13  C23  C33  )  (13-25) 


(13-26) 


(13-27) 


(13-28) 
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From  Equation  (13-28): 


cfc, 

pT  p 

C1  '*'2 

1 

CO 

O 

H  — 

u 

1  0 

0 

cjc, 

cjc2 

cjc3 

= 

0  1 

0 

cjc, 

cjc2 

cjc3_ 

0  0 

1 

(13-29) 


The  matrix  C  is  called  normal  when  CTC  produces  unity  along  principal  diagonal  and 
orthogonal  when  the  off-diagonal  elements  are  0.  When  both  conditions  occur,  the  matrix  is  called 

ortho-normal. 


Equating  the  diagonal  elements: 

CjCj  =  1 ;  CjC2  =  1 ;  CjC3  =  1 


(13-30) 


Using  Equation  (13-24)  and  Equation  (13-25): 


CjC, 


(Cn  C21  c31) 


"11 

“1 

-21 

-31 


=  C?,  +  cl,  +  Cf,  =  1 


^2^2  “  (C ,2  C22  C32)  - 


-12 

-22 

-32 


^-12  +  ^-22  ^-32  ~~  ^ 


C3C3  ~  (C13  C23  ^-33) 


-13 

-23 

1 

-33 


—  Cj3  +  c23  +  c33  —  1 


(13-31) 


(13-32) 


(13-33) 
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As  a  result  of  round-off  errors  in  computations,  the  C  matrix  becomes  increasingly  non-ortho- 
gonal.  If  we  assume  that  these  vectors  are  normalized  by  the  method  described  previously  in 
Section  12,  we  have: 

|Cj|  =  |C2|  =  |C3|  =  1  (13-34) 

where  referring  to  (13-34),  we  constructed  a  set  of  new  vectors,  C] ,  C2,  and  C3  such  that 
=  iCu  +  jC21  +  kC31 
C2  =  *C12  +  jC22  +  kC32 
C3  =  iC13  +  jC23  +  kC33 
and  where,  based  on  the  above  equations: 

icj =  ycfi  +  c21  +  c31 
ic2i =  yc;2  +  c22  +  c32 
1^3! =  yc?3 + c23  +  c33 

As  shown  in  Figure  13-3,  Cj  and  C2  are  not  perpendicular  to  each  other,  at  angle  0.  If  we 
want  to  make  CT  and  C2  to  be  perpendicular  to  each  other,  we  will  need  to  orthogonalize 
Cj  and  C2. 


Figure  13-3.  Non-Orthogonal  Set  of  Column  Vectors 
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Draw  a  perpendicular  line  from  Q  (the  tip  of  C2)  to  Cj,  and  call  the  intersection  P. 

By  definition  of  the  dot  product: 

C2C1  =  C2  Cl  =  ICjIICjI  cos0  (13-35) 

The  vector  OP  is  the  component  of  C2  along  the  direction  of  Cj.  Denoting  the  unit  vector 

Ci 

along  Cj  by  r  ,  we  have: 

_  ,  ,  C, 

OP  =  |C2|  cos0 

[Cil 

"  1  zl  cos  |C,I  |C,| 

=  |C2|  |C,|  cose 

I'-'ll 

T  C, 

=  Cl  C2  — ^  (13-36) 

ic/. 

using  (13-35). 

Denote  the  vector  PQ  by  C 2. 

Then,  from  Figure  13-3,  by  a  vector  addition, 

C2  =  OP  +  C2  (13-37) 

or 

C2  =  C2  -  OP 
=  C2  -  (Cjc,) 

by  using  Equation  (13-36).  (13-38) 

Therefore,  C2  computed  by  Equation  (13-38)  using  the  known  values  of  C2  and  Cj  is  perpendic¬ 
ular  to  Cj. 
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We  started  with  Cj,  C2,  and  C3,  and  made  C'i  perpendicular  to  Cj.  The  next  step  is  to  find 
C3  (replacing  C3),  and  make  C3  perpendicular  to  both  Cj  and  C2. 

Define  Cj  by  C]  =  aCj  noting  that  C/ is  parallel  to  C1;  and  likewise  define  C2  by 
C2  =  bC2  noting  that  C2  is  parallel  to  C2,  where  a  and  b  are  scalar  constants.  (13-39) 

Referring  to  Figure  13-4,  and  based  on  geometry: 

OP  =  Ci  +  C2 

=  aC!  +  bC2  (13-40) 

Also  referring  to  Figure  13-4: 

C3  =  OP  +  PQ 


=  (Cj  +  C2)  +  PQ 

=  aCj  +  bC2  +  C3  (13-41) 
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It  follows  from  Equation  (13-41): 

C3  =  C3  -  aCj  -  bC2 

Now  C3  is  perpendicular  to  C2. 

That  is, 

C2  •  C3  =  0  or  (Cafe's  =  0. 

Using  Equation  (13-42)  in  Equation  (13-43): 

(ci)TCj  =  (Ci)T(C3  -  aC,  -  bci) 

=  (C2)TC3  -  a(ci)TC,  -  b(ci)TC2 

=  (C2)TC3  -  b(ci)TC2 
=  0 

since  (C2)  Cj  =  0  because  we  previously  made  C2  perpendicular  to  Cl5  so  that 

C2  •  Cj  =  (ci)TC!  =  0. 

It  follows  from  Equation  (13  -44) : 

b(ci)Tci  =  (ci)Tc3 


(13-42) 

(13-43) 


(13-44) 

(13-45) 

(13-46) 


Now,  we  want  to  make  C3  perpendicular  to  Cj  (as  well  as  perpendicular  to  C2)  or  we  require: 


Cj  •  C3  =  0  and  C}  C3  =  0 


(13-47) 
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Using  Equation  (13-42): 

Cj  C3  =  C}  (c3  -  aC,  -  bC2) 

=  Cj  C3  -  aCj  Cj  -  bCj  C2 

=  0  (13-48) 


Since  Cj  is  perpendicular  to  C2 ,  Cj  C2  =  0,  we  have  from  Equation  (13-48): 


a  C[  Cj  =  Cj  C3 


(13-49) 


or 


C[c3 

a  =  (13-50) 


Substituting  Equation  (13-46)  and  Equation  (13-50)  into  Equation  (13-42): 
C3  =  C3  -  aCj  -  bC2 


=  C 


3 


ClC3p 

CTC,1-1 


(c2)tc3 

- f - 

(C2)  C2 


(13-51) 


Equation  (13-51)  may  be  evaluated  in  terms  of  known  values  of  Ci,  C3,  and  C2  given  by  Equa¬ 
tion  (13-38). 


Next  we  show  that  C3  i.  C2,  using  Equation  (13-51)  for  C3: 

CTC 


(c2)  c3  =  (c2) 


r  -  1  r  -  r' 

3  C[C,  1  {C2fc2  2 


(13-52) 
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It  follows: 


/  ,XT  ,  ,  ,vT  Cj  C,  ,  ,.T 

(c2)  c3  =  (c2)  c3  -  ^  (c2)  c, 


-  (ci)T  C3  -  (c2)t  c3 


0 


(c;)T  c, ,  ,vT  , 
7-TT-1  N  C2 

(C’2)  C2 


(13-53) 


T 

because  02-LCj,  and  therefore  (C2)  CT  =  0. 


Equation  (13-53)  shows  that  C2  is  perpendicular  to  C3. 

It  follows  that  Cj,  C2,  and  C3  (computed  from  given  Cj,  C2,  and  C3)  are  mutually  orthogonal. 

The  new  orthogonal  set  Ci,  C2,  andC3  (denoting  Ci  by  Ci),  would  be  slightly  off-normal 
because  of  the  way  it  was  derived.  Thus  we  may  have  to  re-normalize,  which  makes  it  slightly  non- 
orthogonal.  So,  the  process  may  have  to  be  repeated  until  both  normality  and  orthogonality  criteria 
are  met  within  error  tolerance  requirements.  Usually  only  a  few  cycles  will  be  sufficient. 
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SECTION  14 

THEOREM  OF  CORIOLIS 


Consider  two  frames,  Frames  A  andB,  with  a  common  origin  at  0.  Frame  B  rotates  with  respect 
to  Frame  A  with  angular  velocity  wab-  A  vector  Rop  =  R  may  be  expressed  in  Frame  A  by 


RA  =  IX  +  JY  +  KZ 


(14-1) 


and  in  Frame  B  by 

Rb  =  ix  +  jy  +  kz. 


(14-2) 


In  the  following  equation,  we  use  a  shorthand  notation  P( )  for  ^  ( ),  and  Pa  for  ^  ( )  with  the 

differential  increment  observed  in  Frame  A,  and  Pb  for  ^( )  with  the  differential  increment 
observed  in  Frame  B,  and  so  on. 


j 

Now,  taking  “j~RA  relative  to  Frame  A,  and  denoting  it  by  PaRa: 

PaRa  =  I#  +  +  K^  +  X§  + 

A  dt  dt  dt  dt  dt  dt 


_  TdX  ,  TdY  ,  x^dZ 

-Idr  +  JdT  +  KdT 


(14-3) 


since  unit  vectors  I,  J,  and  K  are  fixed  in  Frame  A  and  do  not  vary  with  time  in  Frame  A.  Using 
Equation  (14-2), 


;dx  ,  .dy  ^  udz  ,  „ di  ,  dj  ^  „dk 


P  pB  _  jOx  +  srl  +  kSZ  +  v§  +  v^  + 

*ak  *dt  +  +  x~t  +  y^  + 


dt 


dz 

dt 


dt 


dt 


dt 


(14-4) 


In  this  case,  the  unit  vectors  i,  j,  and  k,  which  are  fixed  in  Frame  B,  rotate  relative  to  Frame  A, 
and  therefore  are  time-varying  if  viewed  from  Frame  A. 

Now,  consider  di  in  Equation  (14-4). 
dt 

Referring  to  Figure  14-1,  since  vector  i  is  a  unit  vector,  it  cannot  change  its  magnitude. 
However,  it  can  change  its  direction  relative  to  Frame  A  because  i  is  fixed  in  Frame  B,  which  rotates 
relative  to  Frame  A.  For  an  infinitesimal  angular  displacement,  the  tip  of  the  i  vector  moves  on  the 
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plane  parallel  to  the  j  -  k  plane.  We  can  see  this  in  the  following  way.  If  the  frame  is  rotated 
infinitesimally  about  the  i-axis,  the  direction  of  the  i-vector  is  not  changed. 

If  the  frame  is  rotated  infinitesimally  about  the  z-axis,  the  tip  of  the  i-vector  moves  (to  the  left) 
to  the  direction  parallel  to  the  j-axis.  If  the  frame  is  rotated  infinitesimally  about  the  j-axis,  the  tip 
of  the  i-vector  moves  (downward)  to  the  direction  anti-parallel  to  the  k-direction.  Therefore,  we 
may  decompose  Ai  caused  by  the  angular  displacement  of  i  of  Frame  B  relative  to  Frame  A  in  terms 
of  its  displacements  in  the  j  and  k  directions. 

Referring  to  Figure  14-1,  Frame  A  (with  unit  vectors  I,  J,  K)  coincides  initially  with  Frame  B 

d0z 

(with  unit  vectors  i,  j,  k).  We  see  that  the  angular  velocity  wz  =  — about  the  Z-axis  causes 

A0Z  d0y 

Ai  =  jA0z  =  j  At  =  jwzAt  during  ^t,  and  the  angular  velocity  wy  =  in  Y-axis 

A0y 

causes  Ai  =  —  k  0y  =  —  k  At  =  -  k  wy  At  during  At-  Summing  the  two  components, 
we  have: 

Ai  =  jwz  At  —  kwy  At  (14-5) 

Dividing  by  At  and  taking  the  limits,  we  have 

H  =  jwz-kwy  (14-6) 
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Figure  14-1.  Incremental  Rotation  of  Frame 


The  right  side  of  Equation  (14-6)  is  also  equal  to  wB  x  i  as  shown  below,  where  the 
superscript  B  in  wB  indicates  that  the  components  of  w  are  expressed  in  Frame  B.  Since 

wB  =  (wX;  wy>  wz)T  and  i  =  (1,  0,  0)T  , 

we  have 


j  k 

wy  wz  =jWz-kwy  (14-7) 

0  0 

using  the  definition  of  the  vector  cross-product. 

It  follows  that: 

^  =  wB  x  i  =  jwz  -  kwy  (14-8) 
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Similarly,  we  can  show  that: 

f  =  WB  x  j  (14-9) 

|  =  wBxk  (14-10) 


Now,  referring  to  the  last  three  terms  of  the  right  side  of  Equation  (14-4),  and  using 
Equations  (14-8),  (14-9),  and  (14-10): 

xf  +  yf  +  Zf  =  x(»B  x  i)  +  y(wB  X  j)  +  z(WB  x  k) 

=  wB  x  ix  +  wB  x  jy  +  wB  x  kz  (since  x,  y,  z  are 

scalars) 


=  wB  x  (ix  +  jy  +  kz) 

=  wB  x  Rb  (14-11) 


Note  wB  is  the  angular  velocity  of  Frame  B  relative  to  Frame  A  with  components  expressed 
in  Frame  B). 

Now, 


PbRB  =  PB(ix  +  jy  +  kz) 


=  jdx  +  jdy  kdz 
1  dt  +Jdt  +  kdt 


(14-12) 


because  i,j,  k  are  fixed  in  Frame  B,  and  therefore  ^  =  -  0  and  =  0.  Note  that  the  right 

dt  dt  dt 

side  of  Equation  (14-12)  is  equal  to  the  first  three  terms  of  the  right  side  of  Equation  (14-4). 

So,  substituting  Equation  (14-12)  and  Equation  (14-11)  into  Equation  (14-4): 

PaRb  =  PbRb  +  wBB  x  Rb  (14-13) 


14-4 


NSWCDD/MP-99/17 


In  Equation  (14-13),  the  vector  RB  may  be  expressed,  after  appropriate  transformations,  in  ei¬ 
ther  Frame  A  or  Frame  B,  (see  Section  15).  So,  dropping  the  superscript  B  in  Equation  (14-13), 

PaR  =  PbR  +  wab  x  R  (14-14) 

which  is  called  the  Equation  of  Coriolis  or  the  Theorem  of  Coriolis. 

What  the  Equation  of  Coriolis  implies  is  that  the  differentiation  (with  respect  to  time)  of  a  vec¬ 
tor  in  one  frame  is  not  equal  to  the  differentiation  of  the  same  vector  in  another  frame  that  is  rotating 
with  respect  to  the  first  frame.  And  to  obtain  the  value  of  PAR  in  terms  of  PgR,  we  must  add  a  correc¬ 
tion  term  wAB  x  R  (which  incorporates  the  effect  of  rotation  of  Frame  B  relative  to  Frame  A)  to 
PgR  as  shown  in  Equation  (14-14). 

Readers  who  are  interested  in  the  geometrical  approach  in  the  derivation  of  Equation  (14-14) 
may  find  it  in  other  text  books  such  as  “ Mechanics ”  by  Keith  R.  Symon  (published  by  Addison 
Wesley).  Some  may  find  the  geometrical  approach  difficult  to  follow,  while  others  may  not.  For 
this  reason,  an  analytical  approach  is  presented  here  to  assist  the  comprehension  in  view  of  the  im¬ 
portance  of  the  theorem. 
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SECTION  15 


MATRIX  FORM  OF  THE  THEOREM  OF  CORIOLIS 


The  Theorem  of  Coriolis  in  vector  form  (which  we  derived  previously  in  Section  14)  in  which 
Frame  B  in  rotating  relative  to  Frame  A  with  the  angular  velocity  w^b  is  given  by: 


par  =  pbr  +  wAB  x  R. 


(15-1) 


where  PAR  means  d_i>  observed  in  Frame  A, 
dt 

— R 

and  PBR  means  dt  observed  in  Frame  B. 


Expressing  components  of  Equation  (15-1)  in  Frame  B: 

[par]b  =  [pbr]b  +  w®B  X  rb 


(15-2) 


Now,  referring  to  the  left  side  of  Equation  (15-2)  (using  the  notations  previously  explained  in 
Section  5), 

[paR]“  =  cj  [paR]a 


=  cB  [pra 


=  CA  [P  [CBRB] 


since  RA  =  c£  RB 


CB 


[p  cj]  rb  +  cA  [prb] 


(15-3) 


remembering  that  P( )  is  the  operator  for  ^  ( ),  and  recalling  from  calculus  that  the  chain  rule  given 
below  as: 


dt 


[xy]  = 


d.x 

dtX 


,  dy 

y  +  dT 


which  is  valid  in  the  operation  of  Equation  (15-3)  as  well. 
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It  follows  from  Equation  (15-3): 


[Par]B  -  Cl  [pc( 


rb  +  cB  C£  prb 


rb 


+  prb 


(15-4) 


since  CB  C£  =  Cg  =  I. 

Returning  to  Equation  (15-2)  and  recalling  that  the  vector  wBB  XR  may  be  expressed  in 
matrix  form  (as  explained  in  Section  1)  by: 

wbb  x  RBo[wB|]RB  (15-5) 


T 

where  wBB  =  (wxwywy)  and 


wBK  = 
WAB 


0 

wz 

-Wv 


-Wz 

0 

Wv 


-wx 

0 


Substituting  Equation  (15-4)  and  Equation  (15-5)  into  Equation  (15-2), 
CB  [peg]  Rb  +  PRb  =  PRb  +  [wB|J  Rb 

noting  that  (PBR)B  =  (PRB)B  =  PRB. 

It  follows  that: 


Since  Equation  (15-7)  has  to  be  valid  for  all  RB,  we  conclude: 
C?  P[CJ]  =  wB§ 


(15-6) 


(15-7) 


(15-8) 


Now,  pre-multiplying  (from  the  left)  both  sides  of  Equation  (15-8)  by  (CB)_1  =  (CB)T  =  CB 

[because  CB  is  a  Coordinate  Transformation  Matrix  (between  two  orthogonal  frames),  and  thus  is 
an  Orthogonal  Matrix],  we  have: 


C$  CB  P 


CA 


=  CB  WAB 


,BK 


(15-9) 
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Since  Cg  C®  =  =  I ,  from  Equation  (15-9): 

PC*  =  C*  w“ 

5  cb  =  C*  w*f  (15-10) 


Equation  (15-10)  is  the  matrix  form  of  the  Equation  of  Coriolis. 
Since 

ppA _  d  cA  _.CA(n+l)-CA(n) 

~  dt  Ub  At 


(15-11) 


where  Cg  (n)  denotes  the  value  of  Cg  at  n  At,  where  At  is  the  computation  interval,  we  have  from 
Equations  (15-10)  and  (15-11): 

[n  +  1]  =  Cg[n]  +  Cg[n]  ^[n]  At  (15-12) 

Equation  (15-12)  suggests  that  we  can  update  the  orientation  of  the  Eye  frame  (e.g.,  Frame  B) 
relative  to  Head  frame  (e.g..  Frame  A)  if  we  know  the  angular  velocity  wAB  which  drives  the  eye 
ball. 


wab  §*ves  t^ie  angul31  velocity  of  Frame  B  (e.g.,  Eye  frame)  relative  to  Frame  A  (e.g..  Head 
frame)  with  the  components  expressed  in  Frame  B  in  a  skew-symmetric  form  of  the  3  x  3  matrix 
given  in  Equation  (15-5)  and  repeated  below: 


w 


bk  _ 

AB  ~ 


-wz 

0 

wx 


wY 

-Wx 

0 


(15-13) 


The  vector  equivalent  of  Equation  (15-13)  is 


w 


B  _ 
AB  ~ 


WX 

Wy 

wz 


(15-14) 
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The  angular  velocity  of  B  relative  to  A  with  components  expressed  in  A-frame,  wgB ,  may 
be  obtained  by  a  simple  transformation: 

-iA 

wx 

Wy  =  Cj 
wz 

r  iB 

wx 

=  Cg  WY  (15-15) 

wz 

Then,  the  Skew  Symmetric  Matrix  wgg  may  be  constructed  from  wgB  obtained  by  Equa¬ 
tion  (15-15)  using  Equation  (15-13). 

Next,  we  want  to  derive  the  transformation  equation  to  find  w®^  (with  components  resolved 
in  Frame  B)  from  wg^  (with  the  components  resolved  in  Frame  A),  and  conversely  to  find  wgB 
from  w®£. 

From  Equation  (15-10) 

i  CB  -  C*  =  c§  (15-16) 

Since  A  and  B  are  entirely  arbitrary,  exchange  A  and  B  in  Equation  (15-16): 

C®  =  C®  wgK  (15-17) 

Taking  the  transpose  of  both  sides  of  Equation  (15-17): 

[cb]T  =  [c»  w  «]T 

T  T 

=  [wgj]  [c®]  (see  Appendix  A) 
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Now,  using  Equation  (15-13): 


0  ~wz  WY 
wz  0  ~wx 
-wy  wx  0 

■  0  wz  -Wy 
-wz  0  wx 
wY  ~wx  0 


0 

wz 

-Wy 


=  —  W 


AK 

BA 


-W7  Wv 

0  "WX 

wx  0 


Also,  we  know: 


Using  Equation  (15-19)  and  Equation  (15-20)  in  Equation  (15-18): 


Now  for  the  vector  angular  velocity: 
wab  =  -  WBA 


It  follows  that: 


wAK  = 
WAB 


-  WAK 


BA 


Substituting  Equation  (15-23)  into  Equation  (15-21): 


c§  =  [wab]  C 


(15-19) 

(15-20) 

(15-21) 

(15-22) 

(15-23) 

(15-24) 


Since  A  and  B  are  entirely  arbitrary,  exchanging  A  and  B  in  Equation  (15-24): 
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Now, 

C»C*  =  I 

Since  A)  =  0  ,  differentiating  Equation  (15-26)  using  the  chain  rule: 
£  [c J  Cg]  =  CB  CA  +  d  C* 

=  0 


Substituting  Equation  (15-25)  and  Equation  (15-24)  into  Equation  (15-27): 


WBK  pB  p A  .  pti  AK  PA  _  A 
WBA  UA  ^B  +  UA  WAB  MB  ~  U 


-^B  „7AK  pA  _ 


Using  Equation  (15-26)  in  Equation  (15-28): 

w  BK  _  pB  w  AK  pA 
WBA  UA  WAB  UB 


=  cB 


W 


AK 

AB 


CA 


Using  Equation  (15-23)  in  Equation  (15-29): 


w 


BK  __  pB  ^  AK  f-'A 


BA 


BA  B 


which  transforms  t0  w  ba  • 

Since  A  and  B  are  entirely  arbitrary,  exchange  A  and  B  in  Equation  (15-30): 


WAK  _  pA  wjjk  pi 
WAB  UB  W  AB  CA 


,BK  r’B 


which  transforms  w®|  to  w^g 


(15-26) 


(15-27) 


(15-28) 


(15-29) 


(15-30) 


(15-31) 
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SECTION  16 
QUATERNIONS 

Quaternion  has  been  an  unmixed  evil  to  those  who 
have  touched  them  in  any  way. 

Lord  Kelvin  (1892) 


Lord  Kelvin’s  remark  notwithstanding,  quaternions,  as  far  as  their  applications  to  eye  move¬ 
ments  are  concerned,  are  useful. 

The  use  of  quaternions  greatly  facilitates  the  derivation  of  equations  for  the  rotation  angle  and 
rotation  vector  in  the  Euler’s  Theorem  of  total  equivalent  rotation  and  the  derivations  of  Rodrigues 
equation  as  indicated  in  the  latter  half  of  this  report.  Otherwise,  the  derivations  of  these  equations 
seem  to  be  almost  too  complicated  to  be  tractable  if  we  use  conventional,  non-quaternion  algebra. 

Historically,  the  quaternion  is  the  result  of  the  search  for  a  “Three-Dimensional  Complex 
Number.”  A  complex  number  z  =  x + iy  can  represent  a  vector  r  in  the  plane.  The  complex  numbers 
furnished  the  algebra  for  vectors. 

But  complex  numbers  are  applicable  only  when  all  vectors  lie  on  the  same  plane.  To  treat  vec¬ 
tors  in  three-dimensional  space,  an  analogue  of  complex  numbers  in  three-dimensional  space  be¬ 
came  necessary.  The  mathematicians  in  the  first  half  of  the  nineteenth  century  searched  for  the 
three-dimensional  complex  numbers  and  associated  algebra.  This  search  led  to  the  invention  of 
quaternions  by  William  Rowan  Hamilton  in  1843,  which  inspired  the  emergence  of  vector  algebra 
in  the  latter  part  of  the  nineteenth  century.  It  should  be  noted  that,  historically,  quaternion  algebra 
preceded  vector  algebra,  not  the  other  way  around. 

16.1  QUATERNION  ALGEBRA 

A  quaternion  q  is  a  number  of  the  form  of 

q  =  q0  +  iqi  +  jq2  +kq3  (16-1) 

For  example, 

q  =  2  +  3i  +  5j  +  6k  (16-2) 

in  which  i,  j,  k  play  roles  somewhat  similar  to  i  in  the  complex  the  number  z=a+ib.  The  real  part 
qo  of  the  quaternion,  q,  is  called  scalar  part  s(q),  and  (iqi  +  jq2  +  kq3)  is  called  vector  part  v(q). 
Thus  Equation  (16-1)  is  sometimes  written  as  q=s(q)  +  v(q). 
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The  conjugate  of  q  denoted  by  q*  is  obtained  by  making  the  vector  part  of  q  negative.  That 
is, 


q*  =  qo  -  (iqi  +  jq2  +  kq3) 

=  qo  -  iqi  -  jq2  kq3  (16-3) 

The  three  coefficients  of  the  vector  part  are  rectangular  Cartesian  coordinates  of  a  point,  while 
i,  j,  k  are  unit  vectors  along  the  three  orthogonal  axes. 

The  unit  vectors  i,  j,  k  obey  Hamilton’s  rules  which  are  as  follows  : 

ii  =  i2  =  jj  =  j2  =  kk  =  k2  =  -1  (16-4) 

(Note  the  similarity  of  the  preceding  rule  to  that  of  the  conventional  complex  number,  where 
i  =  /-  1  has  the  result  i2  =  -1.) 


ij  =  k 

ji  =  -k 

jk  =  i 

kj  =  -i 

ki  =  j 

ik  =  -j 

(16-5) 

Note  that  the  sign  convention  for  Equation  (16-5)  is  somewhat  similar  to  that  of  vector  cross- 
products,  where,  for  example,  i  x  j  =  k  and  j  x  i  =  -k,  etc. 

Two  quaternions  p  and  q  may  be  added  or  subtracted  in  a  way  similar  to  that  of  complex  num¬ 
bers.  For: 


p  =  po  +  ipi  +  jp2  +  kp3 
q  =  q0  +  iqi  +  jq2  +  kq3 

We  can  find  p+q  and  p-q  simply: 

p  +  q  =  (p0  +  qo)  +  i(pi  +  qi)  +  j(p2  +  q2>  +  k(p3  +  q3) 
p  -  q  =  (po  -  qo)  +  Kpi  -  qi)  +  j(P2  -  q2)  +  k(p3  -  q3)  (16-6) 

Next,  we  want  to  find  out  if  pq  is  equal  to  qp,  that  is,  if  the  commutative  law  of  multiplication 
holds.  Now: 

pq  =  (po  +  ipi  +  jp2  +  kp3)  (q0  +  iqi  +  jq2  +  kq3)  (16-7) 
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If  we  carry  out  the  operations  of  Equation  (16-7)  using  the  rules  of  Equation  (16-4)  and  Equa¬ 
tion  (16-5),  we  get: 

pq  =  (poqo  -  Piqi  -  P2q2  -  P3q3) 

+  i(poqi  +  piqo  +  P2q3  -  P3q2) 

+  j(poq2  -  piq3  +  P2qo  +  p3qi) 

+  k(p0q3  +  Piq2  -  P2qi  +  P3qo)  (16-8) 

while 

qp  =  (q0  +  iqi  +  jq2  +  kq3)  (p0  +  ipi  +  jp2  +  kq3) 

=  (qoPO  -  qipi  -  q2P2  -  q3P3) 

+  Kqopi  +  qiPo  +  q2P3  -  q3P2) 

+  j(qoP2  -  qiP3  +  q2Po  +  q3Pi) 

+  k(q0p3  +  qlP2  -  q2pi  +  q3p0)  (16-9) 

Compare  Equation  (16-8)  and  Equation  (16-9).  Although  real  (scalar)  parts  are  equal,  not  all 
the  components  of  vector  parts  (i,  j,  k  parts)  have  the  same  signs.  We  see  that: 

pq#qp  (16-10) 

That  is,  the  commutative  law  of  multiplication  does  not  hold  for  the  quaternion  multiplications,  un¬ 
like  the  complex  number  multiplications. 

16.2  QUATERNION  OPERATION  ON  VECTORS 

Quaternions  may  be  used  to  rotate  a  vector  as  well  as  to  change  the  length  of  a  vector.  Consider 
a  vector  r  =  ix  +  jy  +  kz  .  Now  assume  r  is  operated  on  by  a  quaternion  q=qo  +  iqi  +  jq2  +  kq3 
from  the  right  to  become  r'  =  ix’  +  jy’  +  hz’.  That  is, 

q^r'  (16-11) 

or 

(qo  +  iqi  +  jq2  +  kq3)  (ix  +  jy  +  kz)  =  ix'  +  jy'  +  kz'  (16-12) 
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After  carrying  out  the  operation  of  Equation  (16-12),  using  Hamilton’s  Rule,  we  get 
-(xqi  -yq3  -  zq3)  +  i(xq0  +  zq2  -  yq3) 

+  j(yqo  -  zqi  +  xq3) 

+  k(zq0  +  yqi  -  xq2) 

=  0  +  ix'  +  jy'  +  kz'  (16-13) 

Equating  the  coefficients  of  Equation  (16-13)  for  real  parts  and  i,  j  and  k  parts,  we  get: 


xqi  -  yq2  -  zq3  =  0 

(16-14) 

xq0  +  zq2  -  yq3  =  x' 

(16-15) 

yq0  -zqi  +  xq3  =  / 

(16-16) 

zqo  -  yqi  -  xq2  =  z' 

(16-17) 

We  now  have  four  equations  for  four  unknowns  qo,  qi,  q2,  and  q3.  These  are  the  known  values 
of  the  initial  coordinates  x,  y,  z  of  r,  and  the  final  coordinates  x',  y',  z'  of  r',  which  is  sufficient  to 
solve  for  q0,  qi,  q2  and  q3. 

Thus  if  we  want  to  rotate  and  elongate  a  given  vector  r  (x,  y,  z)  into  a  new  vector  r'(  x',  y',  z'), 
we  can  do  so  mathematically,  by  the  operation  qr  =  r'  given  by  Equation  (16-12)  after  first  solving 
for  qo,  qi,  q2,  q3  using  Equation  (16-14). 

16.3  LISTING’S  LAW  IN  TERMS  OF  QUATERNION 

The  use  of  quaternion  leads  to  a  very  simple  formulation  of  Listing’s  Law  (see  Section  3). 

A  quaternion  q  may  be  written  as: 

q  =  q0  +  iqi  +  jq2  +  kq3.  (16-18) 

The  vector  part  (iqi+jq2  +  kq3)  of  q  points  in  the  direction  of  the  axis  of  eye  rotation.  Listing’s 
Law  says  that  the  axis  must  lie  in  the  j-k  plane  (called  equatorial  plane).  This  simply  means  that 
the  coordinate  qi  in  i-direction  is  zero. 

Thus,  we  can  now  state  the  effect  of  Listing’s  Law  very  simply:  The  torsion  of  the  eye  in  any 
position  is  determined  by  a  quaternion  whose  first  vector  component  qi  is  zero,  or  by 

q  =  qo  +  jq2  +  kq3  (16-19) 
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16.4  EXAMPLE  OF  QUATERNION  ALGEBRA  USING  HAMILTON’S  RULES 

For  two  quaternions,  v  and  q,  denoted  by 

v  =  0  +  ix  +  jy  +  kz  (with  the  scalar  part  equal  to  0) 

q  =  qo  +  iqi  +  jq2  +  kq3 , 

vq  =  (ix  +  jy  +  kz)  (q0  +  iqi  +  jq2  +  kq3) 

=  ixq0  +  i2xqi  +  ijxq2  +  ikxq3 
+  jyqo  +  j»yqi  +  j2yq2  +  jkyq3 

+  kzqo  +  kizqi  +  kjzq2  +  k2zq3 
=  ixq0  -  xqi  +  kxq2  -  jxq3 
+  jyq0  -  kyqi  -  yq2  +  iyq3 
+  kzq0  +  jzqi  -  izq2  -  zq3 

It  follows: 

vq  =  -  (xqi  +  yq2  +  zq3) 

+  i(xq0  +  yq3  -  zq2) 

-  j(xq3  -  yqo  -  zqi) 

+  k(xq2  -  yqi  +  zq0) 

Thus,  vq  generates  a  new  quaternion  denoted  by  p  =  po  +  ipi  +  jp2  +  kp3 
where 

PO  =  -  (xqi  +  yq2  +  zq3) 

Pi  =  (xqo  +  yq3  -  zq2) 

P2  =  -  (xq3  -  yqo  -  zqi) 

P3  =  (xq2  -  yqi  +  zq0) 


(16-20) 


(16-21) 


(16-22) 


(16-23) 
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SECTION  17 

EULER’S  ROTATION  VECTOR  AND  RELATIONSHIP  BETWEEN  QUATERNIONS 

AND  ROTATION  MATRICES 


As  we  have  shown  in  Section  16,  a  vector  r  operated  on  by  q  from  the  left  resulting  in  qr 
which  in  turn  generates  another  vector  r'  so  that  qr  =  r'. 

Similarly,  a  vector  v  =  ix+jy  +  kz  pre-multiplied  by  q-1  and  post-multiplied  by  q  generates 
an  new  vector,  v'  =  ix'  +  iy'  +  kz'. 

That  is, 


v'  =  q"1  vq  (17-1) 

It  follows: 

(ix'  +  jy'  +  kz') 

=  (Qo  -  i(3i  “  J^2  ~  kq3)(ix  +  JY  +  kz)(q0  +  *qi  +  jq2  +  kq3)  (17-2) 
In  the  rotation  matrix  notations: 

v'  =  Cv  (17-3) 


or 


'C 

C 

c 


11 

Cl2 

Cjs 

21 

C22 

C-23 

31 

C32 

^33 

IYI 


(17-4) 


where  the  matrix  C  here  is  used  in  the  sense  of  Section  8,  rather  than  Section  5. 

The  above  equation  is  equivalent  to: 

ix'  =  i(Cn  x  +  C12y  +  C13) 

j  y'  =  j  (C21  X  +  C22  y  +  c23) 

kz'  =  k  (C31  x  +  C32  y  +  C33)  (17-5) 
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If  we  carry  out  the  algebra  of  Equation  (17-2)  by  using  Hamilton’s  Rule  (discussed  in 
Section  16.1)  we  get  using  shorthand  notations  Qjj: 

v'  =  ix'  +  j  y'  +  kz' 

=  i(Qn  x  +  Q12y  +  Q13z) 

+  j  (Q21 x  +  Q22  y +  Q23  z) 

+  k  (Q31  x  +  Q32  y  +  Q33  z)  (17-6) 

in  which: 

Q11  =  q?  +  q?  -  qi  -  qi  =  cn 

Ql2  =  2(q0q3  +  qj  q2)  =  C12 

Q13  =  2  (qi  q3 -  qoq2)  =  ^13 

Q21  =  2  (qi  q2  —  qo  *13)  =  ^21 

Q22  =  q?  -  q?  +  qi  -  q?  =  c22 

Q23  =  2  (qo  qi  +  q2  q3)  =  ^23 

Q31  =  2  (qo  q2  +  qi  Q3)  =  ^31 

Q32  =  2  (q2  q3 ~  q0qi)  =  c32 

Q33  =  q?  -  qi  -  qi +  qi  =  c32  (17-7) 

Note  that  Equation  (17-6)  may  also  be  written  as: 


Y' 

Qn 

Ql2 

Ql3 

X 

y' 

z' 

— 

Q21 

Q22 

Q23 

y 

(17-8) 

Q31 

Q32 

Q33 

J 

z 
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Comparing  Equation  (17-8)  with  Equation  (17-4): 


Cll 

C12 

C13] 

Q11 

Ql2 

Ql3] 

^21 

C22 

^23 

= 

Q21 

Q22 

Q23 

C31 

C32 

C33 

Q31 

Q32 

Q33 

where  Qy  is  given  in  Equation  (17-7)  in  terms  of  the  elements  qo,  qi,  q2,  and  q3  of  a  quaternion  q. 
From  Equation  (17-7)  and  Equation  (17-9): 

cn  =  Qn  =  q?  +  q?  -  q2  ~  q? 

c22  =  Q22  =  q§  -  q?  +  q2  -  q? 

C33  =  Q33  =  q?  -  q?  -  ql  +  qi  (17-10> 

It  follows  from  Equation  (17-10): 

C11  +  C22  +  C33  =  3  ql  -  q\  -  -  qj 

=  3  ql  +  q§  -  (qg  +  qj  +  qj  +  qj) 

=  4qJ-  1  (17-11) 

because  q^  +  q  j  +  q2  +  q|  =  1,  (as  derived  in  Equation  (19-10)).  It  immediately  follows  from 
Equation  (17-11): 

qo  =  2  +  +  ^22  +  ^33)2 


(17-12) 
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Also,  from  Equation  (17-7)  and  Equation  (17-9): 


^-"23  ^32  —  ^23  Q32 

=  4qiq0 

(17-13) 

_  ^-23  ~  ^32 

Ql  4q0 

(17-14) 

^•31  -  ^13  =  Q31  “  Q 13 

=  4q2q0 

(17-15) 

_  c3,  -  c13 

q2  4  q0 

(17-16) 

^•12  ~  ^21  =  Ql2  -  Q21 

= 4  q0  q3 

(17-17) 

_  -  c21 

q3  “  in. 

(17-18) 

Thus,  we  now  have  expressions  for  qo,  qi ,  q2,  and  q3  in  terms  of  the  elements  Qj  of  the  rotation 
matrix  with  Equation  (17-12)  in  which  qo  is  given  in  terms  of  Cn,  C22,  and  C33. 

Conversely,  if  want  to  express  Qj  in  terms  of  qo,  qi,  q2,  and  q3,  they  are  given  in 
Equation  (17-7). 


Using  Equation  (17-12)  through  Equation  (17-18),  the  vector  part  v(q)  of  the  quaternion  q  may 
be  expressed  in  terms  of  the  elements  Qj  of  the  rotation  matrix  as  follows: 


v(q)  =  i 


C23  C32 


C3I  C13 


(17-19) 


where  qo  is  given  by  Equation  (17-12): 

q0  =  2  t1  +  Cl1  +  C22  +  C33)2 
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Now,  we  want  to  express  the  rotation  vector  along  the  Euler  axis  in  terms  of  the  unit  vector  u 
and  its  components  ux,  uy,  and  uz  along  X,  Y,  and  Z  axes.  From  Section  19,  on  “The  Angle  of  Quater¬ 
nion,”  we  have: 


u,  =  ux  = 


u2  =  %  = 


u3  =  uz  =  (17-20) 

sin  j 

Using  Equation  (17-14),  Equation  (17-16),  and  Equation  (17-18)  in  Equation  (17-20): 


u7  = 


_  1 

^23  C32 

.  <}) 
sm  j 

2(l  +  Cn  +  C22  +  C33)- 

_  1 

^31  -  C13 

.  4> 
Sm2 

l 

2(1  +  cn  +  c22  +  c33)2 

_  1 

C12  -  ^21 

smf 

2(1  +  Cn  +  C22  +  C33)2 

(17-21) 


where  y  is  the  one  half  of  the  rotation  angle  (J)  about  the  Euler  axis. 


below: 


Now,  we  want  to  express  Equation  (17-21)  in  terms  of  <[>,  not  y,  utilizing  the  identity  shown 


J  (l  +  C11  +  C22  +  c33)2  =  sin(t> 


which  will  be  derived  at  the  end  of  this  section,  and  given  in  Equation  (17-34). 


(17-22) 
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Using  Equation  (17-22)  in  Equation  (17-21): 


ux  = 


_  ^23  ^32 


2  sintj) 


uy  = 


^31  C13 

2  sine}) 


uz  = 


^12  ^21 
2  sin  4> 


(17-23) 


Equation  (17-23)  shows  the  components  of  the  unit  vector  u  =  iux + juy  +  kuz  along  the  Euler 
axis  in  terms  of  Cy  of  the  rotation  matrix  and  Euler’s  principal  angle  <{>.  Although  it  is  relatively 
easy  to  derive  here  by  means  of  quaternion  algebra,  it  would  have  been  extremely  complicated,  had 
we  not  used  quaternions. 

Next,  we  derive  Equation  (17-22). 

From  a  standard  mathematical  table: 

cos  2  a  =  2  cos2  a  -  1  (17-24) 

sin  2  a  =  2  sin  a  cos  a  (17-25) 

From  Equation  (17-11): 

Cjj  +  C22  +  C33  +1  =  4  q  q  (17-26) 


Since  q0  =  cos^  by  definition  (see  Section  19),  from  Equation  (17-26): 

C11  +  c22  +  c33  +  1  =  4  cos2  J  (17-27) 

<t> 

From  Equation  (17-24),  replacing  2a  by  <J>  so  that  a  =  — : 
cos  ({)  =  2  cosz  -y  —  1 

fy  (t) 

cos({)  +  1  =  2  cos  ^ 

2(cos<|>+l)  =  4  cos2  |  (17-28) 
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Using  Equation  (17-28)  in  Equation  (17-27): 

Cjj  +  C22  +  C33  +  1  =2  costj)  +  2 


(17-29) 


This  gives  an  important  equation: 

Cjj  +  C22  +  C33  =  2  cos<J)  +  1  (17-30) 

which  says  that  the  trace  (sum  of  the  diagonal  elements  of  the  rotation  matrix)  is  equal  to  twice  the 
cosine  of  the  rotation  angle  about  the  the  Euler  axis  plus  1. 

Referring  to  the  right  side  of  Equation  (17-29),  and  using  the  first  equation  of  Equation  (17-28) 
for  cos  4> : 


(17-31) 


ignoring  the  negative  root. 

Using  Equation  (17-31)  in  Equation  (17-29): 

(C„  +  C22  +  C33  +  if  =  (2  cos 4>  +  2)2  =  2  cos^  (17-32) 
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From  Equation  (17-25)  replacing  2  a  by  (f)  so  that  a  =  y : 

.  (})  6 

sin  4>  =  2  sin  j  cos  j  (17-33) 

Substituting  Equation  (17-32)  into  the  left  side  of  Equation  (17-22),  and  using  Equa¬ 
tion  (17-33)  we  have: 

=  sin  4>  (17-34) 


(l  +  Cn  +  C22  +  C33)2  =  ^sin  2  ^cos 


which  is  what  we  wanted  to  show. 
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SECTION  18 

EULER’S  ROTATION  VECTOR  RATE 


As  shown  in  Section  2,  the  Euler’s  Rotation  vector  u  (u  for  unit  vector)  has  the  same  direction 
cosines  both  in  the  Eye  Frame  and  the  Head  Frame.  That  is. 


uE  =  uH 


It  is  also  true  that 


uE  =  CE  uH 


It  follows: 


CE  uH  =  uH 


Differentiating  Equation  (18-3)  with  respect  to  time: 

dCH  h  ,  rE  duH  _  duH 
dt  U  +  Lh  dt  "  “d T 


(18-1) 


(18-2) 


(18-3) 


(18-4) 


or 


uH  =  (I  -  Cg)  4*5 


dt 


Referring  to  Section  15: 


dt 


A  pA  —  pA  WB* 
^  ’-'B  WAB 


where  we  used  *  in  place  of  k  in  superscript. 
Identifying  H  with  A  and  E  with  B, 

|cg  =  Cgwgg 

Taking  the  transpose  of  Equation  (18-7): 


-  (c? 

=  [whe]  [c: 


E* 

HE 


Tr  ilT 


(18-5) 


(18-6) 


(18-7) 


(18-8) 
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T 

Now,  [Cg]  =  Cg  for  Orthogonal  Matrix.  Therefore, 


ArH 

dt^E 


A  rE 

dt  Ch 


and,  referring  to  Section  1  and  other  sections,  we  have 


F* 

WHE  = 


0 

W7 


-W7  WY ' 

0  -wx 


-WY  wx  0 


Therefore, 


T 

'  0 

wz 

— w  y" 

F* 

[wheJ 

-wz 

0 

wx 

wY 

— wx 

0 

0  -wz  wY  ‘ 
wz  o  — wx 
-wY  wx  o 


=  —  w 


E* 


HE 


A  matrix  A  is  called  skew-symmetric  if  A  =  -AT  such  as  Equation  (18-11). 

Substituting  Equation  (18-9)  and  Equation  (18-11)  into  Equation  (18-8): 

ApE  _  _  WE*  pE 

WHE 

Substituting  Equation  (18-12)  into  Equation  (18-5): 

duH 


=  (i -eg) 


dt 


Substituting  Equation  (18-3)  into  Equation  (18-13): 

duH 


-W^UH  =  (I-Cg)^ 


(18-9) 


(18-10) 


(18-11) 


(18-12) 


(18-13) 


(18-14) 
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As  an  intermediate  step  to  find  ^  in  terms  of  and  u,  we  want  to  find  a  3  x  3  matrix  M  in  terms 
of  0  and  u  with  the  following  characteristics: 


1) 

MtM  =  I 

(18-15) 

2) 

|M|  =  1 

(18-16) 

3) 

Mu  =  u 

(18-17) 

4) 

Mn  +  ^22  +  ^33  =  trace  M 

=  1+2  cos  4> 

(18-18) 

Equation  (18-15)  requires  that  M  be  an  Orthogonal  Matrix.  Equation  (18-16)  follows  from 
Equation  (18-15).  Equation  (18-17)  and  Equation  (18-18)  complete  the  requirements  for  M  to  be 
identical  to  the  Rotation  Matrix  C. 

The  following  Equation  (18-19)  meets  the  above  requirements  from  Equation  (18-15)  through 
Equation  (18-18): 

(See,  for  example,  Space  Craft  Attitude  Dynamics  by  Peter  C.  Hughes,  John  Wiley,  1986) 

M  ^  coscj)  I  +  (1  —  coscj))  uuT  -  sin<j)u  *  (18-19) 


where 


T 

uu  = 


UX 

uY 

UZ 


ux  UY  uz 


ux  ux  Ux  Uy  ux  uz 

UY  ux  UY  UY  UY  UZ 
UZ  UX  UZ  UY  UZ  UZ 


(18-20) 


and  u*  is  the  skew-symmetric  form  of  the  vector  u  =  [ux  uy  uz]t>  similar  to  Equation  (18-10)  for 
w*  in  form,  and  given  by: 


u  = 


0 

uz 

-uY 


-uz  uY 

0  -uX 

ux  0 


(18-21) 
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Similar  to  Equation  (18-11),  it  is  obvious  that 


Expanding  Equation  (18-19): 


m12 

Mia 

1 

0 

0 

m2, 

m22 

m23 

=  coscj) 

0 

1 

0 

M31 

m32 

M33 

0 

0 

1 

+  (1  -  COS(J)) 


Ux  ux  uy  ux  uz 

Uy  ux  Uy  UY  UZ 
uz  UX  uz  uY  Uz 


-  sin  4> 


0 

uz 

-uY 


-uZ 

0 

ux 


u  Y 

-Ux 

0 


(18-22) 


Expanding  the  right-hand  side  of  Equation  (18-22),  and  equating  the  element  by  element: 

Mn  =  (1  —  cos(j>)  Uj  +  cost}) 

M22  =  (1  —  cos4>)  n\  +  cos<|> 

M33  =  (1  —  COS(t>)  U3  +  COS  (f) 

M12  =  (1  —  cos(|))  U]U2  -I-  u3sin(|) 

M21  =  (1  —  cos(J))  u2Uj  —  u3sin(j) 

M23  =  (1  —  cost}))  u2u3  -1-  UjSintf) 

M32  =  (1  —  costj))  u3u2  —  UjSin4> 

M31  =  (1  —  cos4>)  u3Uj  +  u2sin4> 

M13  =  (1  -  cos(}))  UjU3  —  u2sin4>  (18-23) 
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0  -UZ  u Y 

uz  o  -ux 

-uY  ux  0 

=  [uY  uz  —  uz  uY  -ux  uz  +  uz  ux  ux  uY  —  uY  ux] 

=  [0  0  0]  (18-27) 
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*  * 

u  u  = 


0 

-UZ 

c 

k: 

_ i 

0 

- 

uZ 

uY 

uz 

0 

-UX 

uz 

0 

- 

ux 

-Uy 

ux 

0 

-Uy 

ux 

0 

-  UY 

Uy 

Ux 

— UZ 

ux 

* 

ux 

yy 

-u| 

“ux 

uZ 

Uy 

ux 

uz 

Uy 

uZ 

-4 

— 

ux 

(18-28) 


Since  for  the  Euler’s  unit  matrix,  u^  +  Uy  +  u|  =  1,  we  may  express  Identity  Matrix  I  by 


I  = 


1  0  0 
0  1  0 
0  0  1 


ul  +  ul  +  ul  0 

X  Y  L  2  i  2  i  2 

0  UX  +  UY  +  UZ 

0  o 


0 

0 

UX  +  UY  +  UZ 


since 


It  follows  from  Equation  (18-20),  Equation  (18-28),  and  Equation  (18-29): 
u  u  =  uu  —  I 

Now,  from  Equation  (18-19): 

Mt  =  cos  4)1  +  (1  —  cos4>)uut  —  sin4>u* 

=  cos4>It  +  (1  -  cos4>)[uut]  -  sin4>[u*] 

=  cos4>I  +  (1  —  cos4>)uut  +  sin4>u* 

I  =  IT,  [uuT]  =  [uT]  uT  =  uuT  and  Ju*j  =  —  u*. 


(18-29) 


(18-30) 


(18-31) 
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Let  us  determine  if  Equation  (18-15)  is  satisfied.  Using  Equations  (18-31)  and  (18-19): 


MtM  =  cos(j)I  +  (1  —  cos({))uuT  -I-  sinc()u*j  [cos(t)I  +  (1  —  cos(j))uuT  —  sin<}Hi* 
=  cos2(j)I  +  coscj)(l  —  coscj))uuT  —  cos  (j)  sin  ())u* 

+  (1  —  COS  (|))  COS  (})UUT  +  (1  —  cos  (|))2u(uTu)uT  —  (1  —  COS  4>)  u(uTu*j 


+  sin  c()  cos  4>u*  +  sin  4>(1  -  cos  (j))(u*u)uT  -  sin2())u*u* 

After  some  algebra,  Equation  (18-32)  reduces  to: 

MtM  =  cos2<j)I  +  sin24>I 


+  [terms  which  cancel  each  other] 
+  [terms  which  equal  to  zero] 


1 


0 

0 


0  O’ 
1  0 
0  1 


(18-32) 


(18-33) 


This  proves  Equation  (18-15),  establishing  that  M  is  an  Orthogonal  Matrix. 

After  Equation  (18-15)  is  established,  it  is  a  simple  matter  to  prove  Equation  (18-16).  Taking 
the  determinant  of  Equation  (18-15),  and  using  Equation  (18-33): 

|MtM|  =  |Mt|  |M| 

=  |M|  |M|  =  |M|2 

=  |I|  =  1  (18-34) 

It  follows  from  |M|2  =  1 , 

|M|  =  ±  1  (18-35) 

Referring  to  Equation  (18-19),  M  is  a  continuous  function  of  <[),  and  evaluating  |M|  at  <j>  =  0  : 

M  =  |I  +  0  -  0|  =  |I|  =  1  (18-36) 
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Therefore,  we  conclude 
|M|  =  +  1 


(18-37) 


which  proves  Equation  (18-16). 


To  prove  Equation  (18-17),  we  use  direct  multiplication  of  M  given  by  Equation  (18-19)  by  u 
and  using  Equation  (18-24)  and  Equation  (18-26): 


Mu 


cos(})I  +  (1  -  cos4>)uuT  -  sine))u* 


u 


=  coscfm  +  u  (uT  u)  -  costjm  (uTu)  -  sin  <j>  (u*u)J 
=  COs4>U  +  u  -  costfrn  -  0 

=  u  (18-38) 

Finally,  to  prove  Equation  (18-18),  we  use  the  first  three  equations  of  Equation  (18-23): 

Trace  M  =  Mu  +  M22  +  M33 

=  (1  —  COS(f))Uj  +  cos<)>  +  (1  —  cos  (J>)  u^  +  COS  (}) 

=  +  (1  —  coscj))  U3  -I-  COS({) 

=  Uj  +  Uj  +  U3  —  COS<()(uj  +  U  2  +  U3)  +  3cos4> 

=  1  —  cost})  +  3  cost}) 

=  1  +  2cos<|>  (18-39) 


since  u \  +  +  u|  =  1  for  unit  vector. 

In  summary,  we  have  confirmed  that  M  is  identical  to  the  Rotation  Matrix  ,  or 

Cy  =  M  =  cos  4>I  +  (1  —  cose)))  uuT  —  sin  (j)u*  (18-40) 
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Substituting  Equation  (18-40)  into  the  right  side  of  Equation  (18-14)  for  c£ : 


*]  u  =  |l  -  uuT  -  cos  (j)  (i  -  uuT)  -  u  *  sin  4>  J  ^ 

=  [i  -  uuT  +  cos  (j)I  —  cos  (j)  uuT  +  u*  sin  <J>j  ^ 

=  ^  “  U(“T^)  +  C0S4>^  ~  cos<|m(uT^)  +  u*  sin  (J)  ■ 


(18-41) 


Now,  u  is  d  unit  vector.  Therefore,  it  cannot  change  its  magnitude,  but  it  can  change  its  direction 
only.  This  means  that  the  infinitesimal  change  represented  by  must  be  at  the  perpendicular  direc¬ 
tion  (right  angle)  from  the  direction  of  u.  That  is, 


du  t  du  ^ 

dT  =  0 


(18-42) 


since  cos  90°  =  0. 


Using  Equation  (18-42)  in  Equation  (18-41): 


“  whe]  u  =  ^  “  cos<J)^  +  u*sin<j>^ 

=  [(1  —  cos(())I  +  u*sincj)]^~j7 

Multiplying  both  sides  of  Equation  (18-44)  by  the  inverse  of  ^(1  —  co: 
the  left: 

^  =  [(1  -  cos  <J))I  +  u*  sin  (j)]  [-  wj^u] 

Now,  using  Equations  (18-10),  (18-11),  and  (18-21),  we  can  show  that: 


(18-43) 


(18-44) 


—  cos<t>)I  +  u  sin  cf>  from 


(18-45) 


E*  *  F 

-  WHEU  =  U  Wgg 


(18-46) 
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or 


'  0 

wz 

—  W  Y 

'ux' 

'  0 

-  uz 

Uy* 

-wx- 

-wz 

0 

wx 

Uy 

= 

uz 

0 

-UX 

W  y 

Wy 

-wx 

0 

uz 

-Uy 

ux 

0 

wz 

'  Wz  Uy  —  WY  UZ' 
-wzux  +  WXUZ 
Wy  UX  —  WX  Uy 


(18-47) 


This  confirms,  for  two  vectors  a  and  b,  axb  =  —  b  x  a  may  be  expressed  a*b  =  —  b*a 
using  a  Skew-symmetric  Matrix. 

Using  Equation  (18-46)  in  Equation  (18-45): 

^  =  [(1  —  cos({))I  +  u  *  sin4>]  1  u*  w^E  (18-48) 


Note  that  w^E  is  a  3  x  1  column  vector,  while  w^E  is  a  3  x  3  mat  rix  equivalent  of  the  vector 
whe  >  which  is  the  angular  velocity  of  the  eye  relative  to  head  with  components  expressed  in  the  Eye 
frame. 


The  inverse  indicated  in  Equation  (18-48)  may  be  expressed  as 
[(1  —  cos(J))I  +  u*sin4>J 


1 

1  -  COS(f> 


I  -  ^  sin  4>u*  +  i  (1  +  coscj>)  (uuT  -  i) 


(18-49) 


by 
I  = 


The  validity  of  the  right  side  of  Equation  (18-49)  may  be  checked  by  multiplying  it  from  the  left 
(1  —  cost])) I  +  u* sin cj>J  and  confirming  that  the  result  reduces  to  the  Identity  Matrix 


T  0  O' 
0  1  0 
0  0  1 


Substituting  Equation  (18-49)  into  Equation  (18-48): 


1 

1  —  COS({) 


<j)  [uV]  -  i  (1  +  COS<|))  u* 


w 


HE 


(18-50) 
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where  we  used  the  relation 

^  (1  +  cos  4>)  (uuT  —  i)  U*  =  i  (1  +  COS(j))  [u(uT  u*)  —  u  * 


=  —  i(l  +  cosc}))  u* 


since  uTu*  =  [0  0  0]  by  Equation  (18-27). 
It  follows  from  Equation  (18-50): 


du 


1  *  1 


*  1  a..  *  .  * 


uu  1  1  1  *  1  •  1  v 

17  =  i - r  xu  ~  i  cos  d>u  —  sin  <pu  u 

dt  1  -  cos  <j>  [2  2  2  Y 


w 


HE 


1 


1  —  cos(J)  [2 


^■(1  —  cos  (J>)  u*  —  sin  4>  u*  u* 


w 


HE 


From  Equation  (18-52): 


du 

dt 


1*  1  cp**  p 

2U  ~2COt2  U  U  WHE 


in  which  we  used  a  trigonometry  identity: 

_  <l>  _  sincj) 

COt  A  1  t 

2  1  -  cos  <p 

Using  Equation  (18-30)  for  u*u*,  Equation  (18-53)  may  also  be  written  as 


du 

dt 


1  *  1  V  /  T  t\ 

1U  ~  1  COS  j  (uu1  -  I) 


w 


HE 


(18-51) 


(18-52) 


(18-53) 


(18-54) 


(18-55) 


Equation  (18-53)  or  Equation  (18-55)  is  the  time  rate  of  the  unit  vector  along  the  Euler’s  single 
equivalent  rotation. 


The  equation  for  ^-j-  as  well  as  for  the  rotation  angular  rate  (scalar)  about  the  Euler’s  princi¬ 

pal  axis  (see  Section  22)  is  difficult  to  solve.  One  reason  is  that  the  vector  u  is  defined  only  after  a 
rotation  has  taken  place;  consequently,  there  is  some  problem  in  setting  an  initial  value  for  u. 
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SECTION  19 

THE  ANGLE  OF  QUATERNION:  THE  ROTATION  ANGLE  OF  EULER’S  THEOREM 


Consider  a  unit  vector  u  along  the  axis  of  a  single  equivalent  rotation  in  the  context  of  the  Euler’s 
Theorem  that  makes  angle  a  with  both  XA  and  XB  axes,  angle  (3  with  both  YA  and  YB  axes  and 
angle  y  with  both  ZA  and  ZB  axes.  Then,  the  components  ui,  U2  and  U3  of  u  along  X,  Y,  and  Z  axes 
of  both  Frames  A  and  B  are  related  to  a,  (3,  and  y  (see  Section  2): 

ui  =  cos  a  along  X  axis 

U2  =  cos  P  along  Y  axis 

U3  =  cos  y  along  Z  axis. 


The  quaternion  q  may  be  expressed  as: 

q  =  qo  +  iqi  +  jq2  +  kq3  =  qo  +  q  (19-1) 

in  which 

qo  is  called  the  scalar  part  of  q_  denoted  by  s(q),  and 
(iqi  +'jq2  +  kq3)  is  the  vector  part  of  q,  denoted  by  v(q). 

So,  we  may  write: 

q  =  s(q)  +  v(q)  (19-2) 


Assume  Frame  B  (e.g.,  moving  Eye  frame)  is  rotated  by  angle  <j)  from  Frame  A  (reference-pri¬ 
mary  frame)  about  the  Euler  axis. 

Euler  four-parameters  qo,  qi,  q2  and  q3  (which  obey  the  quaternion  algebra  explained  in 
Section  16)  are  defined  by: 

q0  =  cosset)  (19-3) 

q1  =  ut  sin^4>  =  cos  a  sin^<|)  (19-4) 

q2  =  u2  sin  Jcp  =  cosp  sin^tj)  (19-5) 

q3  =  u3  sini<j>  =  cosy  sin^  (19-6) 
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Substituting  Equation  (19-3)  through  Equation  (19-6)  into  Equation  (19-1): 

q  =  cos^-cf)  +  icos  a  sin^4>  +  jcos  |3  sin^<})  +  kcos  y  sin^-tj) 

=  cos^t})  +  (icosa  -I-  j cos |3  +  kcosy)  sin^<J> 


Since 


i  cos  a  +  j  cos  P  +  k  cos  y  =  i  Uj  +  j  u2  +  k  u3  =  u 


we  have  from  Equation  (19-7): 

<j>  .  4> 

q  =  cos  y  +  u  sm  — 

in  which  (|>  (not  is  called  the  angle  of  the  quaternion  q  and  u  its  axis. 
Using  Equations  (19-3),  (19-4),  (19-5),  and  (19-6): 


9o  +  9?  +  9?  +  9? 


=  cos2tt  +  (cos2  a  4-  cos2p  -I-  cos2 y)  sin2 — 


(1)  7  (J) 

=  COS  TV  +  Sin Z  77 


=  1 

since  cos2  a  +  cos2p  +  cos2y  =  1  based  on  Equation  (19-8). 
Adding  the  squares  of  Equation  (19-4),  (19-5),  and  (19-6): 


92  +  q2  +  q2  =  (cos2  a  +  cos2|3  +  cos2  y)  sin2  ^ 


•  2<t> 
=  sm2  y 


It  follows: 


•  <t> 
sm2 


=  (9?  +  92  +  9§) 


1/2 


(19-7) 


(19-8) 


(19-9) 


(19-10) 


(19-11) 


(19-12) 
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If  the  Eye  frame  (Frame  B)  is  rotated  from  the  reference-primary  orientation  by  angle  <j>  about 
the  Euler  axis  with  unit  vector  u,  then  we  represent  the  new  orientation  of  the  Eye  frame  by 

4*  •  <t> 

q  =  cos  2  +  usin-j  (19-13) 

Inverse  of  q  :  q-1 

Previously,  we  defined  the  conjugate  of  q*  of  a  quaternion  q,  given  by 


q  = 

qo  +  *qi  +  Jq2  +  kq3 

(19-14) 

to  be 

q  *  = 

qo  ~  (iqi  +  Jq2  +  ^3)  =  %-  ich  -  -  kq3 

(19-15) 

It  follows: 

qq*  = 

(qo +  iqi  +  jq2  +  kqs)(qo  -  »qi  -  jq2  -  ^3) 

(19-16) 

q*q  = 

(qo  -  iqi  -  jq2  -  kqs)(qo  +  *qi  +  jq2  +  kq3) 

(19-17) 

In  both  Equation  (19-16)  and  Equation  (19-17),  if  we  carry  out  the  algebra  using  the  Hamilton’s 
Rule  explained  in  Section  16,  we  get  the  same  result  as  follows: 

q  q  *  =  q  *  q  =  q§ +  q?  +  qi +  =  q? +  (qf  +  +  q|) 

=  cos2  J  +  sin2^  =  1  (19-18) 

using  Equations  (19-3)  and  (19-11). 

By  definition,  the  inverse  denoted  by  q_1  of  q  must  satisfy: 

qq_1  =  q_1q  =  1  (19-19) 

From  Equation  (19-18)  and  Equation  (19-19),  we  conclude  that: 

q-i  —  q  * 

The  inverse  q_1  of  q  is  equal  to  the  conjugate  q*  of  q  if  we  use  the  definition  given  in 
Equations  (19-3),  (19-4),  (19-5),  and  (19-6).  That  is, 

q'1  =  q*  =  qo  -  fai  +  jq2  +  ^3)  (19-20) 
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Now  consider  a  quaternion  operation  on  a  vector  v  such  that  qvq  1  results  in  a  new  vector  v'. 
That  is, 

v  =  qvq-1  (19-21) 

Now,  consider  a  Rotation  Matrix  C  that  performs  a  similar  operation  on  vector  v  to  yield  v'  such 
that 

v'  =  Cv  (19-22) 


t - 

;> 

C11 

C12  C13 

'vr 

V2 

V3 

w  J 

C21 

C31 

^22  C23 

C32  ^33 

> 

v2 

v3 

It  is  shown  in  Equation  (17-12)  of  Section  17, 

q0  =  2  t1  +  Cn  +  C22  +  c33)2 
From  Equation  (19-3)  and  Equation  (19-24): 

cos  i  4>  =  i  (i  +  CU  +  C22  +  C33)2 

Squaring  Equation  (19-25): 

COS2  ±  4>  =  \  (l  +  C11  +  C22  +  C33) 

Using  trigonometry,  identify  from  a  Math  Table: 
n  6  1  +  coscb 

C0S  2  =  — T“ 

and  substituting  Equation  (19-27)  into  Equation  (19-26): 

1+2COS4>  =  \  (1  +  C„  +  C22  +  C33) 


(19-23) 


(19-24) 


(19-25) 


(19-26) 


(19-27) 


(19-28) 


It  follows: 


C11  +  C22  +  C33  —1  +  2  cost}) 


(19-29) 


Thus,  the  sum  of  the  diagonal  elements,  called  “trace,”  of  the  Rotation  Matrix  C  is  equal  to 
(1  +  2  cos  0).  From  Equation  (19-29): 


<})  =  cos 


A  (c„  +  c22  +  C33  -  1) 


(19-30) 
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in  which  4>  is  the  rotation  angle  about  the  single  equivalent  rotation  axis,  which  moves  Frame  A  (pri¬ 
mary-reference  frame)  to  Frame  B  (current  Eye  orientation). 


To  confirm  that  §  given  in  Equation  (19-29)  is  indeed  the  rotatiori  angle  about  the  Euler  axis, 
consider  a  rotation  by  0z  about  Z-axis.  Then,  the  Rotation  Matrix  from  Frame  A  to  Frame  B  is,  as 
shown  in  Section  5  by: 


fCn 

C12 

C13 

cos  0Z 

sin  0Z 

O' 

C-21 

^22 

^23 

= 

-sin  0Z 

cos  0Z 

0 

C31 

^32 

C33 

0 

0 

1 

(19-31) 


Applying  Equation  (19-29): 

Cji  +  C22  +  C33  —  COS  ©2  +  COS  ©2  +  1 

=  1  +  2cos0z 


Comparing  Equation  (19-29)  with  Equation  (19-32): 
1  +  2cos(f)  =  1  +  2cos0z 


or 


—  0Z 


as  expected. 

Similarly,  for  the  rotation  about  X-axis  by  0x: 


fCn 

C12 

C13] 

1 

0 

0 

C21 

C22 

C23 

0 

cos  0X 

sin0x 

C31 

^32 

C33 

a 

0 

-sin0x 

cos0x 

Cn  +  C22  F  C33  —  1  +  cos0x  +  cos0x 
=  1  +  2  cos  0X 

Comparing  Equation  (19-29)  with  Equation  (19-36): 

1  +  2  cos  <(>  =  1  +  2  cos  0 


(19-32) 

(19-33) 

(19-34) 

(19-35) 

(19-36) 

(19-37) 


or 


4>  =  0X 


(19-38) 


as  expected. 
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SECTION  20 
RODRIGUES  VECTOR 


According  to  Euler’s  Theorem,  any  sequence  of  rotations  of  a  rigid  body  (such  as  an  eye-ball 
or  the  eye-coil  frame  fixed  on  it),  that  has  one  point  fixed  at  the  origin  of  the  coordinate  frame  can 
be  achieved  by  a  single  equivalent  rotation  through  some  angle  cp  (called  Euler’s  principal  angle) 
about  some  axis  (called  Euler’s  axis)  that  passes  through  the  same  fixed  point  at  the  origin  of  the 
coordinate  frame. 

Referring  to  Section  2,  suppose  the  unit  vector  u  along  Euler’s  axis  makes  an  angle  a  with  the 
Xh  axis  of  Head  frame  and  the  Xe  axis  of  Eye  frame;  an  angle  3  with  the  Yh  axis  and  Ye  axis;  and 
an  angle  y  with  the  Zh  axis  and  Ze  axis. 

Then  the  components  ux,  uy,  and  uz  of  the  unit  vector  along  the  X-axis,  Y-axis  and  Z-axis 
of  both  frames  are 

ux  =  cos  a  ;  uY  =  cos|3  ;  uz  =  cosy  (20-1) 


and 

u  =  iux  +  juY  +  kuz 

In  Section  19,  we  defined  the  parameters  for  the  quaternion  q  by, 
q  =  q0  +  iqi  +  jq2  +  kq3 

by 

4> 

q0  =  cos j 

.  <i>  .  <i> 

qj  =  cos  a  sin =  uxsm-^- 

(b  (b 

q2  =  cos  p  sin  —  =  uysin— 

.  4>  .  4> 

q3  =  cosysm—  =  uzsm-y 

where  (j)  is  Euler’s  principal  angle  for  the  single  equivalent  rotation. 


(20-2) 


(20-3) 


(20-4) 
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Note  that  Equation  (20-4)  is  defined  in  terms  of  not  (|). 

Rewrite  Equation  (20-3)  as 

q  =  s(q)  +  v(q)  (20-5) 

where 

s(q)  =  (scalor  part  of  q)  =  q0 

v(q)  =  (vector  part  of  q)  =  iq2  +  jq2  +  kq3 


Using  the  definition  of  Equation  (20-5),  we  express  the  vector  part  of  q  as: 

v(q)  =  iqj  +  jq2  +  kq3 

cb  (b  <h 

=  IUX  sin  2  +  jUy  sin  j  +  kuz  sin-^ 

=  (iux  +  juy  +  kuz)siny 

=  u  sm  j 


(20-6) 


Equation  (20-6)  shows  that  the  vector  part  (iqj  +  jq2  +  kq3)  is  equal  to  the  unit  vector  u  along 
<b  (h 

Euler’s  axis  scaled  by  sin  y  •  We  modify  v(q)  by  dividing  it  by  cos  -^  =  q0  and  call  it  the 
Rodrigues  Vector  p  (or  some  call  it  Gibbs  Vector).  The  Rodrigues  Vector  p  is  defined  by: 


*  u  sm^r 

p=— ^ 

cosf 


=  u  tan 


(20-7) 


where  4  implies  “is  equal  to”  by  definition. 

From  Equation  (20-7)  with  Equation  (20-2): 

iuxsm4-  +  iuvsin~-  +  kuzsm4- 

p  =  2  J  y  2  . - l  (20-8) 

COSy 
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It  follows  using  Equation  (20-4): 


.92  ,  i  93 
Jtt  +  Ktt- 
J9o  qo 


or 


ran 

qo 

q2 

qo 

q 3 
qo 


(20-9) 


(20-10) 


Note,  in  Equation  (20-7),  p  goes  to  infinity  when  the  denominator  cos^  =  0  or  ^  =  90° 

or  <j>  =  180°.  This  singularity  at  (j)  =  180°  of  pis  a  limitation  in  the  use  of  the  Rodrigues  vector, 
while  referring  to  Equation  (20-4),  the  magnitude  of  the  quaternion  parameters  q0  qj  q2  and  q3 
cannot  exceed  unity.  An  advantage  of  the  Rodrigues  parameters  over  the  quaternion  parameter  is 
that  the  former  has  three  parameters  while  the  latter  has  four  parameters. 


However,  the  singularity  at  4>  =  180°  poses  no  problem  for  eye  movement  analysis  because 

this  situation  equivalent  to  the  rotation  of  the  primary,  reference  eye  position  180°  to  the  back  of 
the  head,  is  physiologically  impossible  movement. 


In  Section  17,  we  found  that: 

q0  =  i  a  +  c„  +  c22  +  c33)'/2 


qi  = 


c23  c32 


2  (1  +  Cn  +  C22  +  C33)V2 


q2  = 


C31  C13 


2  (1  +  C11  +  C22  +  C33)*/2 


93  - 


C12  c21 


2  (1  +  Cjj  +  C22  + 


(20-11) 

(20-12) 

(20-13) 

(20-14) 


Using  the  above  equations,  we  can  express  the  Rodrigues  parameters  in  terms  of  Cy  of  the 
Rotation  Matrix  C. 
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It  follows  from  Equation  (20-10),  using  Equation  (20-11)  to  Equation  (20-14): 


o  -  qi  - 

^23  ^32  .... 

R32  R23 

(20-15) 

Px  “  q0  " 

i 

+  cn  +  c22  +  c33 

1  +  Rn  +  R22  +  R33 

^31  “  C13  _ 

R13  _  R31 

(20-16) 

Py  "  q0  ~ 

T 

+  cn  +  c22  +  c33 

1  +  Rn  +  R22  +  R33 

Cl2  “  ^21  _ 

R21  -  R12 

(20-17) 

pz  “  q0  “ 

1  +  cn  +  c22  +  c33 

1  +  Rn  +  R22  +  R33 

eliminating  the  square  roots  present  in  the  denominators  of  Equation  (20-11)  to  Equation  (20-14). 
For  the  meanings  of  Ry  vs  Cy,  see  Section  5. 

In  Section  17,  we  found  the  relationship  between  the  elements  Qj  of  the  Rotation  Matrix  C  and 
the  parameters  q0,  q](  q2,  and  q3  of  quaternion  q: 

cll  C12  C13]  fq0  +  qi  "  q2  “  2(clo  ^3  +  qi^)  2(q1q3-q0q2) 

C2]  C22  C23  =  2(qiq2-q0q3)  q20  -  q?  +  q\  -  qj  2(qoqi  +  q2q3)  (20-18) 

c3l  C32  C33J  2(q0q2  +  q1q3)  2(q2q3-q0q1)  qg  -  qj  -  q2  +  q2 

From  Equation  (20-10),  we  have 

qi  =  PiOo  ;  q2  =  p2qo ;  q3  =  p3qo  (20-19) 

(Pi  =  Px>  P2  =  Py’  P3  =  Pz) 

Previously  in  Section  19,  Equation  (19-10),  we  found  that 

9o  +  q?  q2  ql  =  i  (20-20) 

Substituting  Equation  (20-19)  into  Equation  (20-20): 

qo  +  p5qo  +  p2  Po +  P3  q?  =  q§  +  p5  +  p2  +  pf)  =  1  (2°-21) 


or 


_ 1 _ 

1  +  pj  +  p2  +  p2 


(20-22) 
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Using  Equation  (20-19)  and  Equation  (20-22)  in  several  elements  of  Equation  (20-18): 

cn  =  <io  +  qi  -  ol  -  4 

=  9?  +  pjqJ  -  pjil  -  p3q» 

=  qg  (1  +  P?  -  Pl  -  pf)  (20-23) 

c>‘  ■  1  +  p;+1p^p^ «  +  Pi  -  Pi  -  Pi)  <20'24> 

C12  =  2  (qoq3  +  qiq2) 

=  2  (qoP3<lo  +  PiQc^^o) 

=  2  qj  (p3  +  Plp2) 

C'2  =  1  +  P2  +  p2  +  p2  (P3  +  P1P2)  (20‘25) 

Similarly,  the  rest  of  Qj  may  be  found  in  terms  of  P]  p2  and  p3.  The  results  are: 


Cll 

^12 

C13] 

C21 

C22 

C23 

= 

C31 

C32 

C33 

*4 

f1+ 

Pi  “ 

P2-P3 

1  +  P2  +  P2  +  p2 


>3  +  PiP2)  2(Pl  p3  -  p2) 
2(pi  P2  P3)  1-P?  +  P2“P3  2(Pl  +  P2P3> 

2(P2  +  Pi  P3)  2(p2P3-pj)  l-p2-p2  +  p2 


(20-26) 


To  determine  Euler  principal  angle  for  single  equivalent  rotation,  we  use  the  previously 
found  relation: 

Cjj  +  C22  +  C33  =  1  +  2  cost})  (20-27) 

in  which  C 1  j  C22,  and  C33  may  be  computed  in  terms  of  Pl  p2  and  p3  from  Equation  (20-26). 
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SECTION  21 

RODRIGUES  VECTOR  DIFFERENTIAL  EQUATION 


In  this  section,  we  derive  the  differential  equation  for  the  Rodrigues  Vector  driven  by  the  angu¬ 
lar  velocity  of  the  Eye  frame  relative  to  the  head  acting  as  input. 


Referring  to  Section  20,  the  Rodrigues  Vector  p  is  defined  by: 
<j> 

p  =  u  tan  j 


(21-1) 


where  u  is  the  unit  Euler’s  Rotation  vector  and  (j)  is  the  (single  equivalent)  rotation  angle  about  Eul¬ 
er’s  Rotation  axis.  It  follows: 


dp  _  d_ 
dt  dt 


using  the  identity  ^  tan  A  =  sec2  A  ^  where  A  is  an  angle. 
In  Sections  18  and  22,  respectively,  we  derive: 


Substituting  Equation  (21-3)  and  Equation  (21-4)  into  Equation  (21-2): 


(21-2) 

(21-3) 

(21-4) 

(21-5) 
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Referring  to  the  last  term  in  Equation  (21-5), 

wTu  =  uTw  =  wx  ux  -I-  WY  Uy  +  wz  uz  =  w  •  u  =  u  •  w 
which  is  a  scalar. 

Thus,  we  can  write: 

uwTu  =  u(uTw)  =  (uuT)w 

Substituting  Equation  (21-6)  into  Equation  (21-5)  for  u(wTu),  and  factoring  out  w: 


(21-6) 


(21-7) 


where  trigonometry  identities 
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By  definition  of  the  Rodrigues  Vector  p, 


*  4*  * 

p  =  tanyu 


(21-8) 


which  is,  in  expanded  form  (of  skew-symmetric  matrix): 


'  0 

“Pz 

Py' 

4> 

=  tan-£ 

'  0 

“UZ 

uY- 

z=z 

Pz 

0 

“Px 

uz 

0 

— UX 

-pY 

Px 

0 

2 

“UY 

ux 

0 

It  follows  from  Equation  (21-8): 

H-) 


p*p  =  tan- 


4*  * 

=  tan-ju 


9  (b  t 

=  tan^u1!!. 


(21-9) 


Using  Equation  (21-8)  and  Equation  (21-9)  in  Equation  (21-7),  we  have 

t  -  [i +  p* +  ppT]f 


(21-10) 


Equation  (21-10)  may  also  be  expressed  in  vector  form 


^  =  2  [w  +  pxw  -I-  p(p  •  w)] 


(21-11) 


in  which  we  used  (ppT)w  =  p(pTw)  =  p(p  •  w),and  p*w  is  replaced  by  the  vector  cross-product 
pxw. 

Equation  (21-10)  or  Equation  (21-11)  is  the  differential  equation  of  the  Rodrigues  vector,  driv¬ 
en  by  the  angular  velocity  of  the  Eye  frame  relative  to  the  Head  frame  with  components  expressed 

in  Eye  frame,  that  is,  driven  by  wHE  . 

Pre-multiplying  both  sides  of  Equation  (21-10)  from  the  left  by  the  inverse  of  |l  +  p*  +  ppTj 
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we  have. 


w  =  2 


I  +  P  +  PP 


In  Section  26,  we  derive  the  following  equation  for  w: 


w 


1  +  pTp 


T  * 

1  -  P 


1  dp 
dt 


For  Equation  (21-12)  to  be  consistent  with  Equation  (21-13),  we  have  to  show  that 

r  *  Tl-1  [I  -  pi 

['  +  p  +  pp  1  =  TTofo 


To  validate  the  above  expression  we  must  show  that: 

HtS  I1  +  p’  +  ppT]  - 1 


(which  is  to  be  proven). 

Referring  to  the  left  side  of  Equation  (21-15): 

[i  -  p']  [i  +  p-  +  pPT] 

T  .  *  .  nr  *  ** 

=  I  +  p  +  PP  -  P  -  P  P  -  P  PP 


(21-12) 


(21-13) 


(21-14) 


(21-15) 


x  *  *  *  T 

=  I  +  pp1  -  p  p  -  p  pp1 


(21-16) 
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Now 


0  ~Pz  Py"  [Px 

p*p  =  Pz  0  “Px  Py 

“Py  Px  0  Pz 


“Pz  Py  +  Py  Pz  0 
Pz  Px  ~  Px  Pz  =  0 
“Py  Px  +  Px  Py  0 


0  ~Pz  Py'  ’  0  “Pz  Py' 
p*p*  =  Pz  0  “Px  Pz  0  -Px 

[“Py  Px  0  J  [“Py  Px  0 


(21-17) 


-p|  -  py  Py  Px  Pz  Px  " 

=  Px  Py  “Pz  “P  y  Pz  Py. 

Px  Pz  Py  Pz  “Py“Px 


[PxPyPz] 


px  1  0 

[PxPyPz]  Py  0  1 

pz  0  0 


'p2  Px  Py  Px  Pz' 

Py  Px  Py  PyPz 
Pz  Px  Pz  Py  Pz 


Px +  Py  +  Pz  0  ® 

0  Px  +  Py  +  Pz  .  °7  9 

0  0  Px  -*■  Py  +  Pz 


Py  -p|  Px  Py  Px  Pz  ' 
=  Py  Px  -PrPy  PJ  Pz. 
Pz  Px  Pz  Py  “Px“Py 


(21-18) 


(21-19) 
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Comparing  Equation  (21-18)  and  Equation  (21-19),  we  conclude: 

p*p*  =  ppT  -  pTpI  (21-20) 

Substituting  Equation  (21-17)  and  Equation  (21-20)  into  Equation  (21-16): 

[i-p*]  [r  +  p'  +  pp7] 

- 1  +  ppt  -  pV  -  (p*p)pt 


=  I  +  ppT  -  ppT  +  pTpI  -  0 


=  I  +  pTp! 


=  (l  +  PTp)l  (21-21) 

Substituting  Equation  (21-21)  into  Equation  (21-15),  we  have 

[I  -  p*]  [I  +  p*  +  ppT]  _  (1  4-  pTp)I  _  1 

1  +  PTP  1  +  PTP  (21-22) 


QED. 
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SECTION  22 

ROTATION  ANGULAR  RATE  ABOUT  THE  EULER  AXIS 


In  Equation  (17-30),  we  found  that 

Cn  +  C22  +  C33  =  1+2  cos  <j)  (22-1) 

where  <()  is  the  rotation  angle  about  the  Euler  Axis  of  the  Eye  frame  relative  to  the  Head  frame,  and 
Cn,  C22,  and  C33  are  the  diagonal  elements  of  the  rotation  matrix  from  Frame  E  to  Frame  H,  as 
shown  below: 


I-H  =  eg  rE 


cn 

C21 

C3j 


C12 

C22 

^32 


Cj3 

^23 

C33 


XH 

Yh 

ZH 


XE 

Ye 

ZE 


(22-2) 

(22-3) 


Note  that  the  trace  of  Cg  is  the  same  as  the  trace  of  C|,  being  the  sum  of  the  diagonal  elements 
that  are  common  to  both  matrices. 

From  Equation  (22-1): 

cost})  =  ^(CH  +  C22  +  C33  -  1)  (22-4) 

By  definition,  Cn  +  C22  +  C33  is  called  the  trace  of  C,  denoted  by  tr  C.  Using  this  notation 
in  Equation  (22-4): 

costj)  =  i[tr  C|  -  l]  (22-5) 


dd) 

We  want  to  find  -37-. 

dt 


Differentiating  both  sides  of  Equation  (22-5)  with  respect  to  time: 

d 


•  1 
““O' dT  =  2 


1- 


d  pE 

dt^H 


(22-6) 
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In  Section  15,  we  found: 


A  pH 
dt  E 


=  cH  wE* 

UE  WHE 


(22-7) 


where  w^g  is  the  matrix  form  of  the  angular  velocity  w|g  of  Frame  E  relative  to  Frame  H,  with  the 
superscript  E  indicating  that  its  components  are  expressed  in  the  Frame  E. 


Taking  the  transpose  of  both  sides  of  Equation  (22-7): 


Now: 

[c?]T  =  eg 

-|T 

0  -wz  wY 
wz  o  -wx 

— WY  WX  0 
— W  y 

wx 

0 

Wz  Wy 

0  -wx 

wx  o 


0  wz 
“wz  0 

Wy  -WX 

0  - 
-  wz 

-Wy 


=  —  W 


E* 

HE 


(22-8) 


(22-9) 


(22-10) 


Substituting  Equation  (22-9)  and  Equation  (22-10)  into  Equation  (22-8): 


d_ 

dt 


—  w 


E* 

HE 


c 


E 

H 


Substituting  Equation  (22-11)  into  Equation  (22-6): 


S,n<^  dT  ~  ~2 


tr 


w 


E* 

HE 


c 


E 

H 


(22-11) 


(22-12) 
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It  follows: 


^  _ 1 tr[wE*  CE1 

dt  2sin<j)  trr  he  lhJ 

T 

Recalling  that  C|  =  [Cg]  ,  and  referring  to  Equation  (22-3): 


(22-13) 


0  — wz  wY 

-L-  tr«  wz  0  -wx 
2smq>  -Wy  wx  q 


11 

^21 

C31 

12 

C22 

^32 

13 

^23 

c33 

(22-14) 


Remembering  that  the  trace  of  a  matrix  is  the  sum  of  the  diagonal  elements  or  (1,1),  (2,2)  and 
(3,3)  elements,  we  have  after  expanding  and  summing  diagonal  elements: 

d(J)  i  r 

d t  2sirT<j>  ' (_Wz °12  +  Wy Ci3)  +  (wz  C2i  ”  wx  C23)  +  (-wY  C31  +  wx C32)] 

=  +  2^  [' WX  (C23  -  C32)  +  Wy  (C31  -  C13)  +  wz  (C12  -  C21 )]  (22-15) 

In  vector  form,  the  matrix  jwj^gj  becomes 


WHE  =  WY  and  WjJjg  =  [wx  wY  wz] 
wz 

In  Section  17,  the  unit  vector  along  the  Euler  axis  was  found  to  be: 

fuxl  C23  “  C32 

u  =  uY  =  -  -1  -  C31  -  C13 


(22-16) 


Y  2  sin  c()  31  13 

uz  C12  -  C21 


(22-17) 


where  Cy  in  Equation  (22-17)  are  the  elements  of  Cg,  which  is  equal  to  the  transpose  of  Cj^ . 
Using  Equation  (22-16)  and  Equation  (22-17): 

WHE  u  =  2ihT6  [Wx  (C23  -  c32)  +  wY  (c31  -  C13)  +  wz  (C12  -  C21)]  (22-18) 
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Comparing  Equation  (22-15)  and  Equation  (22-18),  we  conclude: 


d<J> 

dt 


u 


T 

Note  that  (w|E)  =  [wxwYwz]  is  in  a  vector  form,  while 


w 


E* 

HE 


0 

w. 


-w, 


-wz 

0 

wx 


Wy 

-wx 

0 


(22-19) 


is  in  matrix  form  of  the  angular  velocity  of  the  Eye  frame  relative  to  the  Head  frame. 

Equation  (22-19)  is  a  compact  form  of  Equation  (22-15),  which  expresses  the  Euler  angle  rate 
Equation  (22-15)  expresses  ^  in  terms  of  4>,  and  components  of  w  and  C,  while  Equa- 
dd> 

tion  (22-19)  expresses  -j^-in  terms  of  vectors  w  and  u. 
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SECTION  23 

UNIFIED  FORM  OF  EULER’S  ROTATION  VECTOR  RATE  AND  ROTATION 
ANGULAR  RATE  ABOUT  THE  EULER  AXIS 


In  Section  18,  we  derived  the  following  equation  for  the  Euler’s  Rotation  Vector  Rate: 


where 


u 


du  _  1 
dt  2 


u*  —  cos(^l(uuT 


-I) 


w 


HE 


(23-1) 


u  =  unit  vector  along  the  Euler’s  Rotation  axis. 


u  = 


'  0 

uz 

-uY 


-u. 


u-l 


0  .  -u 
ux  0 


X 


(23-2) 


=  a  Skew-symmetric  Matrix  corresponding  to  vector  u  = 


UX 

uY 

uz 


<})  =  scalar  rotation  angle  about  the  Euler’s  Rotation  Axis. 

WHE  =  angular  velocity  vector  of  the  Eye  frame  (Rotating  frame)  relative  to 
the  Head  frame  (Reference  frame)  with  the  components  expressed  in  the  superscript  Frame  E. 

Also  in  Section  22,  we  derived  the  equation  for  the  Rotation  Angular  Rate  about  the  Euler  Axis 
given  below: 


i  =  7jf  =  (whe)T« 


(23-3) 


in  which  4>  and  (j>  are  scalar  variables. 

To  achieve  the  unified  form,  combining  the  equations  for  u  and  (J),  we  define  a  new  vector 
variable  tj>  by 


<f>  =  4>u 


(23-4) 
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in  which  <{>  has  the  same  direction  as  u  with  its  magnitude  equal  to  $  since  u  is  a  unit  vector.  That 
is,  |<t>|  =  (j). 


It  follows  from  Equation  (23-4): 

<t> 

11  =  x 

<t> 


(23-5) 


Differentiating  Equation  (23-5): 

du  _  1  ,  1  ,  , 

*  -  ^ - 


(23-6) 


Solving  Equation  (23-6)  for  (j)  and  using  Equation  (23-5): 


=  <t>f +  u* 


(23-7) 


By  direct  substitutions  for  u*  from  Equation  (23-2),  it  is  easy  to  show  that 


T  T  *  * 

UU1  “  I  =  u  u 


(23-8) 


or 


T  T  l  *  * 

UU1  =  I  +  U  U 


(23-9) 


where  uT  =  ux  uY  uz] . 


1  • 

Now,  substituting  Equation  (23-1)  for  ^  and  Equation  (23-3)  for  4>  in  Equation  (23-7),  and 
denoting  w^E  by  w  for  simplicity: 


<|>  =  <l> 


iu*  -  icot^j(uuT  -  i) 


w  +  uuTw 


(23-10) 


where  we  used  wTu  =  uTw  corresponding  to  the  dot  product  w  •  u  =  u  •  w  . 
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From  Equation  (23-5): 


u  = 


ux 

uY 

uz 


4> 

§ 

foz 

<}> 


4>x 

4>y 

<}> 


Substituting  Equation  (23-11)  into  Equation  (23-2): 

$Y* 

<i> 

"^X 

4> 

0 


_  1 

0 

(J)Y 

^z 

0 

-4>x 

-4>y 

<l>x 

0 

1 

<i> 

It  follows: 

^  ^  1  .  3|C  |  DjC 

"  =  <t> 


u  = 


0 


— <t>z 

4> 


0 


4> 

~4>y  4>x 

<t>  <1> 


(23-11) 


(23-12) 


(23-13) 
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Using  Equation  (23-9),  Equation  (23-12)  and  Equation  (23-13)  in  Equation  (23-10): 


4>  =  <t> 


i  COt  (f )  ^2  <>*  <)>*  w+(l  +  u*u*)w 


=  1  4>*  w  -  T  cot  (I)  T5  <t>V  w  +  w  +  4>*  <|>*  w 


.2;  4>2 


=  W  +  i<t)*W+^j 


:s] 


i  =  »  +  i<t»*w+i 


,  <t>  sin  cj)  ]  *  * 

‘-j  Tr^  <t>  4>  w 


(23-14) 


where  we  used  a  trigonometry  identity: 


2  1  —  cost}) 


(23-15) 


Equation  (23-14)  is  the  unified,  combined  form  of  u  and  <j)  given  in  Equation  (23-1)  and  Equa¬ 
tion  (23-3). 

We  may  change  Equation  (23-14),  which  is  a  matrix  differential  equation  to  a  vector  differen¬ 
tial  equation,  by  replacing,  in  Equation  (23-14),  (j)*w  by  <j)  X  w  and  <}>*4>*w  by  <{>  x  <|>  x  w  as 
shown  below: 


<i,  =  w  +  I<t)Xw+^(l-f  X  <0  x  w.  (23-16) 

The  above  equation  is  originally  developed  by  John  E.  Bortz  (see  his  article  listed  in  the  Bibliog¬ 
raphy  7).  The  Equation  (23-14)  is  derived  here  by  using  entirely  different  approaches  from  those 
originally  used  in  the  Bibliography  7  to  facilitate  the  comprehensions  by  readers,  consistent  with 
the  developments  in  this  report. 
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SECTION  24 

QUATERNION  DIFFERENTIAL  EQUATION 


Referring  to  Section  19,  we  recall  a  quaternion 


q  =  qo  +  *qi  +  Jq2  +  = 


qo 

qi 

q2 

qs 


=  q0  +  q 


where 


q 


iqi  +  jq2  +  kq3  = 


qi 

q2 

q2 


Equation  (24-1)  may  be  expressed  by  a  column  vector  by:  (see  Section  19) 


4> 

cos  2" 

qo 

qi 

ux  smy 

q2 

•  <l> 

q2 

UySmTj- 

•  <l> 

uzsinj 

qo 

q 


4> 

cos2 


•  <j> 
usin2 


(24-1) 


(24-2) 


where  <J)  is  the  angle  of  quaternion  q,  and  u  is  the  unit  vector  along  the  Euler  Axis  of  single  equivalent 
rotation  and  u  =  iux  +  juy  +  kuz . 
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From  Equation  (24-2),  we  have: 


q  =  u  sin 


4> 


(24-3) 


Differentiating  Equation  (24-3)  with  respect  to  time: 

dq  i  <t>  d(J>  (j)  du 
^  cos  —  u  —  +  sin  -  — 


dt 


From  Section  22, 


dt 


2  dt 


(24-4) 


dcj>  _ 
d T  " 

And,  from  Section  18, 

du 


E  1^ 

wheJ  u 


(24-5) 


±u*-±cot(f)(uuT-l) 


w 


HE 


(24-6) 


where  wj|E  is  the  angular  velocity  of  the  Eye  frame  relative  to  the  Head  frame  with  the  components 
resolved  into  the  Eye  frame. 

Substituting  Equation  (24-5)  for  ^  and  Equation  (24-6)  for  into  Equation  (24-4): 


^  =  {["'he]T»}  §  (cos$)  “ 


+  sin 


(?) 


1*1* 
2  U  “a001 


'( | )  (u  "T  “  *) 


w 


HE 


(24-7) 


It  follows  from  Equation  (24-7),  since 

cot  (I)  = cos  (t)/ sin  (I)  that: 

§  =  isin|u*W  +  icos|(wTu)u 


1  (b  t  i  (1) 

-  2  C0S  2  uu  w  +  2  cos  w 


=  i  sin  ^  u*w  +  i  cos  ^  [(wTu)  u  -  u  (uTw)  +  w 


(24-8) 
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where  we  used  w  for  w|E  for  simplicity. 

Now  uTw  =  wTu  corresponds  to  a  dot  product  u  •  w  =  w  •  u,  which  is  a  scalar  (non- vec¬ 
tor)  quality.  We  can  easily  see  that 

uuTw  =  u(uTw)  =  u(wTu)  =  (wTu)u  (24-9) 

Using  Equation  (24-9)  in  Equation  (24-8): 

do  1  .  <t)  *  i  (b 

dt  =  2  Sm  2  U  "  +  2  C0S  2  W 

We  recall  the  skew-symmetric  form  w*  corresponding  to  the  vector  w  is  given  by: 


(24-10) 


* 

w  = 


o 

wz 

-Wv 


— wz 
0 
Wv 


Wy 

— wx 
0 


(24-11) 


Similarly,  u*  corresponding  to  u  is  given  by: 


* 

u  = 


o 

uz 

-Uv 


0  -Ux 
Ux  o 


(24-12) 


that 


where 


By  direct  substitutions  of  Equation  (24-11)  and  Equation  (24-12),  it  may  be  easily  confirmed 


*  * 

u  w  =  —  w  u 


(24-13) 


Wx' 

Ux 

w  = 

wy 

and  u  = 

Uy 

wz 

Uz 
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Using  Equation  (24-13)  and  Equation  (24-2)  in  Equation  (24-10): 


do  1  .  cb  *  | 

dT  =  "  2  sm  2  W  u  +  2  q0w 


1  *  l 

~  2W  Sm 


q0w 


Now 


\  w*q  +  \  wq0 


dq0 

dt 


d  4> 

dlC0S  2 


1  .  <t>  d4> 
"2Sin2  dT 


Using  Equation  (24-5)  in  Equation  (24-15): 


dq0  1  .  4*  t 

_  =  _  _  sm  -  w?u 


1  .  <J)  T 
-  2  Sin  y  u'w 


if.  <i>  V 

=  -2  “n2u  W 


It  follows  from  Equation  (24-16)  using  Equation  (24-2): 


dqo  1  T 

—  =  —  —  n 1 


dt 


qw 


1  T 
=  —  4  w1' 


w"q 


(24-14) 


(24-15) 


(24-16) 


(24-17) 
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Combining  Equation  (24-17)  and  Equation  (24-14),  and  expressing  in  a  matrix  form: 

dq0 
dq  =  dt 
dt  dq 
dt 


(24-18) 


where 


(24-19) 


24-5 
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Expanding  Equation  (24-19): 

^90  =  ^(°q0  +  -wx9i  “  wyq2  -  wzq3) 

=  ^(wxq0  +  Oqj  +  W2q2  -  wyq3) 

;|q2  =  ^(wyqo  -  w2qi  +  °q2  +  wxq3) 

^q3  =  ^(wzq0  +  wyqi  -  wxq2  +  °q3) 

Rearranging  Equation  (24-20) 

^q0  =  qiwx  -  q2wy  -  q3wz) 

^qi  =  |(q0wx  -  q3wy  +  q2wz) 

^q2  =  f(q3wx  +  q0wy  -  qiwz) 

^q3  =  \[~  q2wx  +  qjwy  +  q0wz) 
Expressing  Equation  (24-21)  in  a  matrix  form: 


dq  d 

qo 

qi 

-qi  -q2  -q3 
qo  -q3  q2 

wx 

Wv 

dt  dt 

q2 

qs 

q3  qo  -qi 
-q2  qi  qo 

J 

Wz 

(24-20) 


(24-21) 


(24-22) 
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SECTION  25 

ANGULAR  VELOCITY  OF  EYE  FRAME  RELATIVE  TO  THE  HEAD  FRAME  IN 

TERMS  OF  QUATERNION  ELEMENTS 

In  this  section,  we  derive  the  angular  velocity  of  the  Eye  Frame  relative  to  the  Head  Frame, 
,  in  terms  of  the  quaternion  elements  qo,  qi,  q2  and  q3. 

Reviewing  some  of  the  notations  we  used  in  the  previous  sections, 


w 


HE 


=  W  = 


wx 

wY 

wz 


[whe]T 

II 

H 

II 

Wx  Wy 

wz] 

'  0 

-wz 

Wy" 

F* 

WHE 

= 

* 

w  = 

wz 

0 

-wx 

— Wy 

wx 

0 

r  nT 

'  0 

wz 

-Wy* 

WHE 

= 

-wz 

0 

wx 

L  J 

Wy 

— wx 

0 

It  follows: 


Kf  - 


—  w 


E* 

HE 


(25-1) 


In  Section  24,  we  derived: 


q0 

0 

-wx 

-wY 

-wz 

4i 

_  1 

wx 

0 

wz 

-Wy 

q2 

2 

Wy 

-wz 

0 

wx 

93 

wz 

Wy 

-wx 

0 

qo 

<Il 

92 

<13 


(25-2) 
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Using  notation  qv  =  <h  in  which  the  subscript  v  refers  to  the  vector  part  of  a  quaternion 

K 

q  =  q0  +  iqj  +  jq2  +  kq3,  we  may  write  Equation  (25-2)  as  shown  in  Equation  (24-18): 


qo  1  0  -wT  [Qo 

qv  ~  2  [w  _W*J  [Qv 


(25-3) 


Expanding  Equation  (25-3): 


2  q0  =  -wTq 


(25-4) 


2  qv  =  wq0  -  w  qv 


(25-5) 


Referring  to  the  term  —  w*qv  in  Equation  (25-5): 


-w*  qv  = 


0  wz  ~wy]  Ri' 
-wz  o  wx  q2 
wy  -wx  o  q3 


wz  q2  -  wy  q3 

wzq,  +wxq3 
wYqi  -  wx  q2 


‘o  -q3  q2l  [wx 

=  q3  o  ~qi  wy 

-<b  qi  o  wz 


(25-6) 


where  qv  is  a  3  x  3  matrix  equivalent  of  a  3  x  1  vector  qv.  See  Section  1.  It  follows  by  putting 
Equation  (25-6)  into  Equation  (25-5)  for  -w*qv : 


2  qv  =  qvw  +  wq0 


(25-7) 


For  any  vectors  A  and  B,  A  •  B  =  B  •  A,orATB  =  BTA.  It  follows  from  Equation  (25-4): 


2  q0  =  -  wTqv  =  -  qTw 


(25-8) 
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Multiplying  both  sides  of  Equation  (25-7)  by  q*  from  the  left: 

2q*  Qv  =  qvqvW  +  q*wq0  (25-9) 

Referring  to  the  term  q^  q*  w  in  Equation  (25-9),  recall  the  vector  identity 

A  x  (B  x  C)  =  B(C  •  A)  —  C(A  •  B) ,  and  that  A*B  is  equivalent  to  the  vector  cross-product 
A  x  B  ,  and  A  •  B  =  B  •  A : 


=  qv(w  •  qv)  -  w(qv  •  qv) 
=  qv(qvw)  -  w(qjqv) 


Now: 


qjqv  —  [919293] 


■qr 

92 

93 


=  9?  +  92  +  93 


=  (9o  +  9?  +  9|  +  9?)  "  9o 
=  1  -  9o 


where  9o  +  9?  +  92  +  9?  =  1  (as  derived  in  Section  19)  is  used. 
Using  Equation  (25-11)  in  Equation  (25-10): 

9v9*w  =  qvqvW  -  w(l  -  qg) 

=  qvqvw  -  (1  -  9§) w 

since  (l  —  q^  is  a  scalar. 

Using  Equation  (25-12)  in  Equation  (25-9): 

2q*  qv  =  qvqvW  -  (l  -  q§)  w  +  q*wq0 


(25-10) 


(25-11) 


(25-12) 


(25-13) 
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Now  pre-multiplying  Equation  (25-4)  by  qv  on  both  sides: 

2qv  Qo  =  ~  qvwTqv  =  -qvqvW  (25-14) 

since  wTqv  =  q3w  for  the  dot  product. 

Putting  Equation  (25-14)  into  Equation  (25-13)  for  qvqvW,  and  transposing: 

2q*qv  +  2qv  %  =  -  (l  -  qj)  w  +  q*wq0  (25-15) 

Now  multiplying  Equation  (25-7)  by  qo  (scalar): 

2q04v  =  q*wq0  4-  wq§  (25-16) 

It  follows: 

q*wq0  =  2q0  qv  -  wqj  (25-17) 

Putting  Equation  (25-17)  into  Equation  (25-15)  for  q*wqo  : 

2q*  qv  +  2qv  q0  =  -  w  +  qgw  +  2q0  qv  -  wqg 

=  -  w  +  2q0qv  (25-18) 


Solving  for  w: 


w 


HE 


=  w  =  2 


qo 


q0dtqv 


(25-19) 


Equation  (25-19)  is  the  desired  results  which  expresses  the  angular  velocity  of  the  Eye  Frame 
relative  to  the  Head  Frame  with  components  expressed  in  the  Eye  frame.  That  is,  w^E  in  terms 


of  the  elements  of  a  quaternion,  qo  and  qv  = 


qi 

92 

93 


and  their  derivative  with  respect  to  time. 
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We  recall  from  Section  17  that 

q0  =  ^(l  +  Cn  +  C22  +  C33) 

=  ^(l  +  Ru  +  R22  +  R33)2 
_  ^23  ~  C32  _  R32  ~  R23 

qi  4qo  4qo 

_  C31  “  C13  _  R13  -  R31 
q2  4q0  "  4q0 

„  _  C12  -  c21  _  r21  -  r12 
q3  4q0  4q0 


(25-20) 

(25-21) 

(25-22) 

(25-23) 
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SECTION  26 

ANGULAR  VELOCITY  OF  EYE  FRAME  RELATIVE  TO  HEAD  FRAME 
IN  TERMS  OF  RODRIGUES  PARAMETERS 


In  this  section,  we  derive  an  equation  to  express  the  angular  velocity  of  the  Eye  frame  relative 
to  the  Head  frame  with  the  components  expressed  in  the  Eye  frame  ( wEE ),  in  terms  of  the  Rodrigues 
parameters  pj,  p3>  and  p3,  which  we  have  discussed  in  Section  20. 


In  Section  25,  we  derived  the  equation  for  in  terms  of  the  quaternion  elements  qo,  qi, 
q2,  and  q3.  The  equation  is  shown  below: 


WHE  =  2 


n 

q°  dt 


dq0  1  ~ 

^  Ov  2qv 


dqv 

dt 


(26-1) 


Here  qv  denotes  the  vector  part  of  a  quaternion  q  =  qo  +  iqi  +  jq2  +  kq3.  That  is 


qv  = 


0i 

q2 

03 


=  »oi  +  j02  +  k03 


(26-2) 


Our  approach  is  to  express  each  term  on  the  right  side  of  Equation  (26-1)  in  terms  of 
Pj,  p2,  and  p3.  Referring  to  Section  20: 


P  = 


Pi 

P2 

_p3 


Si 

q0 
q 2 
Oo 
03 
Oo 


It  follows: 


0i' 

PiOo 

qv  = 

02 

P2q0 

_q3_ 

P3q3 

(26-3) 


(26-4) 
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Using  Equation  (26-4): 

P  •  P  =  Px  +  Py  +  Pz  =  Pi  +  P2  +  P3 

=  q?  +  qj  + 

<lo 


It  follows: 


1  +  p  •  p 


$ 


^5 


+ 


q2  + 


^0 


2 

3 


since  qjj  +  q?  +  +  q3  =  1  (see  Section  19). 

From  Equation  (26-5): 

qo  = - - — T  =  (1  +  P  •  P)_I 

(1  +  P  •  p)2 


Now 


qv  =  iqj  +  jq2  +  kq3 
=  ipiqo +  jp2qo  +  kP3qo 
=  qol'Pl  +  JP2  +  kP3) 

=  q0p 

Substituting  Equation  (26-6)  into  Equation  (26-7): 
qv  =  p(l  +  p  *  p)~J 


Using  the  chain  rule,  the  time-derivative  of  qv  becomes: 


qv=^  =  £[p(1+p‘p)4 


=  ^(1  +  P  *  P)  2  +  P^(l  +  P  •  P)  2 


(26-5) 


(26-6) 


(26-7) 


(26-8) 


(26-9) 
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Referring  to  Equation  (26-9): 

^(1  +  p  '  p)"5  =  -  |  (1  +  p  •  p)“2  (P  •  P  +  P  •  p) 

=  -|(1  +  p  •  p)~’  (2p  •  p) 

=  -  (1  +  p  •  p)~2  (p  •  p)  (26-10) 


Note  that  P  •  p  =  p  •  p  for  dot  product. 

Using  Equation  (26-10)  in  Equation  (26-9): 

qv  =  7u(1  +  P  ‘  p)”'  “  p(1  +  p  ‘  p)_*  (p  '  p) 

Now  from  Equation  (26-6): 

q° "  =  i  (1  +  p ' 

=  -  ^(1  +  p  •  p)”5  (p  •  p  +  p  *  p) 

=  -  (1  +  p  •  p)-^  (  P  *  p) 

Next,  using  Equation  (26-6)  and  Equation  (26-11): 

q04v  =  (1  +  P  •  P)4  +  P  •  P)4-  P(1  +  P  •  P)‘f  (p  • 

=  (1  +  p  •  p)-1  ^  -  p(l  +  p  •  p)“2  (p  •  p) 

Next,  using  Equation  (26-13)  and  Equation  (26-8): 

9oqv  =  -  (i  +  p  •  p)-5  (p  •  p)  p(i  +  p  •  p)  21 


(26-11) 


(26-12) 


(26-13) 


(26-14) 


=  -  (1  +  p  •  p)  2  p(p  •  p) 


(26-15) 
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Finally,  using  Equation  (26-8)  and  Equation  (26-11): 
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-l 


It  follows,  factoring  out  (1  +  p  •  p)  : 


W  =  Wup  =  ,  ,  — - 

HE  1  +  p  •  p 


dp-ox^a 
dt  p  dt 


or  in  matrix  form: 


w  = 


i  * 

1  -  p 


1  +  pTp  L*  ^  J  dt 


dp 


where 


0  — P3  P2 

P3  0  “Pi 

— P2  Pi  0 


p  .  p  =  pTp  =  p2  +  p2  +  p2 


and 


(26-18) 


(26-19) 


(26-20) 

(26-21) 


dp  _  d. 
dt  dt 


TPx" 

Py 

Pz 


(26-22) 
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APPENDIX  A 

DEMONSTRATION  OF  (AB)T  =  (AB)"1  =  BT  AT  =  B"1  A-1 
(FOR  3X3  ORTHOGONAL  MATRICES  A  AND  B) 


A-l/(A-2  blank) 
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DEMONSTRATION  OF  (AB)T  =  (AB)"1  =  BT  AT  =  B"1  A-1 
(FOR  3X3  ORTHOGONAL  MATRICES  A  AND  B) 


This  demonstration  is  straight  forward,  based  on  direct  multiplications. 
For  Rotation  Matrices: 


and 


’  an 

a12 

a13 

'  b„ 

bj2 

bj3 

A  = 

a21 

a22 

a23 

and  B  = 

b21 

b22 

b23 

a31 

a32 

a33 

b31 

b32 

b33 

11 

a21 

a31 

12 

a22 

a32 

13 

a23 

a33 

Bt  =  B"1  = 


bn  b2i  b31 
bi2  b22  b32 


b]3  b23  b 


33 


since  A  and  B  matrices  are  othogonal. 
It  follows: 


(A-l) 


(A-2) 


(A-3) 


all  a12  a13 

bn  b12  b13 

AB  = 

a21  a22  a23 

b2i  b22  b23 

a31  a32  a33 

b3i  b32  b33 

anbn  +  a12b21  +  a13b31 
a2lbll  +  ^22^21  +  ^23^31 
a3ibn  +  a32b21  +  a33b31 


allbl2  a12b22  +  ^3^32 
a2lbl2  +  a22^22  +  ^3^32 
a3lbl2  +  a32b22  +  **33^32 


anbj3  +  a12b23  +  a13b33 
a2lbl3  +  a22^23  +  a23^33 
a3lbl3  +  ^2^23  +  **33^33 


(A-4) 
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It  follows  by  (1)  exchanging  rows  and  columns  of  Equation  (A-4),  (2)  noting  that  the  product 
of  two  Orthogonal  Matrices  is  another  Orthogonal  Matrix,  and  (3)  remembering  that  the  transpose 
of  an  Orthogonal  Matrix  is  equal  to  the  inverse  of  the  original  matrix,  that: 

[AB]t  =  [AB]  ~ 1  = 


allb 

11 

+ 

a12b21 

+ 

a13b31 

a21b 

11 

+ 

a22b21 

+ 

a23b31 

a31bl  1 

+ 

a32b21 

+ 

a33b31 

allb 

12 

+ 

a12b22 

+ 

a13b32 

a21b 

12 

+ 

a22b22 

+ 

a23b32 

a31b  12 

+ 

a32b22 

+ 

a33b32 

allb 

13. 

,+ 

a12b23 

+ 

a13b33 

a21b 

31 

+ 

a22b23 

+ 

a23b33 

a31b 13 

+ 

a32b23 

+ 

a33b33 

Now,  from  Equation  (A-3)  and  Equation  (A-2): 
BtAt  =  B-1A-1 


bll  b21  b31 

all  a21  a31 

= 

b12  b22  b32 

a12  a22  a32 

b,3  b23  b33 

a13  a23  a33 

bnan  +  b21a12  +  b31a13 
bi2an  +  b22a12  +  b32a13 
bi3an  +  b23a12  +  b33a13 


bna2i  +  b21a22  +  b31a23 

b12a21  +  b22a22  +  b32a23 

b13a21  +  b23a22  +  b33a23 


bna3i  +  b21a32  +  b3]a33 
b12a31  b22a32  +  b32a33 

b13a31  +  b23a32  +  b33a33 


(A-6) 


Comparing  Equation  (A-5)  and  Equation  (A-6),  we  conclude: 

[AB]t  =  [AB]-1  =  BtAt  =  B-1A-1  (A-7) 


Equation  (A-7)  is  true  only  for  the  Orthogonal  Matrices. 

For  Non-orthogonal  Square  Matrices,  the  formula  for  the  transpose  and  the  inverse  are  sepa¬ 
rate,  and  Equation  (A-7)  does  not  hold.  That  is, 


[AB]t  =  btat 

(A-8) 

fe 

1 

II 

w 

I 

> 

1 

(A-9) 

But,  generally,  BTAT  ^  B_1A_1. 

By  extension: 

[ABC]1  =  [[AB]C]T  =  Ct[AB]t 

=  CtBtAt 

(A- 10) 
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and 


and 


[ABC]"1  =  [[AB]C]  1  =  C-^AB]”1 
=  C-1  B-1  A-1 

[ABC]1  =  [ABC]"1  =  CtBtAt  =  C-1B-1A-1. 
Now  denote  vector  r  and  rT  by: 


r  = 

X 

y 

and  rT  = 

X  ’ 

y 

z 

z 

-.t 


=  [xyz] 


Then: 


Cr  = 


Cn  C12 

Cl3 

X 

= 

C-21  ^-22 

^23 

y 

C3i  C32 

C33 

z 

‘  Cnx  + 

ci2y  + 

CBz 

= 

C2ix  + 

c22y  + 

C23z 

C3jx  + 

c32y  + 

C33Z 

It  follows: 

[Cr]T  =  CBx  +  C12y  +  CBz  C2jX  +  C22y  +  C22z  C2jX  +  C22y  +  C33Z 
Now, 

,  C11  C21  C31 
rTCT  =  [x  y  z]  |  C12  C22  C32 

C13  C23  C33 


(A-ll) 


(A- 12) 


(A-13) 


(A- 14) 


[CnX  +  C12y  +  CBz  C21x  +  C22y  +  C22z  C31X  +  C32y  +  C33zj  (A-15) 
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Comparing  Equation  (A- 14)  with  Equation  (A- 15),  we  conclude 
[Cr]T  =  rT  CT 

Example:  Application  to  Fick’s  System. 

In  Section  6,  we  found  out  that  for  the  Fick’s  System  of  rotations, 
rE  =  CFj  CFz  CF>  rH 

r  uf2  S,  r 


(A-16) 


(A- 17) 


where 


rE  = 


’  XE  ' 

and  rH  = 

XH 

Ye 

yH 

ZE 

ZH 

Taking  the  transpose  of  Equation  (A-17): 

-  [c£ c*,  ch'  rRlT 


-  C«1 
=  hf  -  [ft'f  [<$  <$■]"] 


[rHf  -  [cj 


-.T 


$[' 


CS  CK 


(A-18) 
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Now: 


-  [XE 

Ye  ze] 

[rf 

=  [XH 

Yh  zh 

-iT 

=  RP. 

hJ 

H 

cf:]t  = R£ 

[cfi]t  =  RfJ 

It  follows  from  Equation  (A- 18)  and  Equation  (A- 19): 
[xe  Ye  ze]  =  [xh  Yh  zh]  Rh  RFj  RF2 
in  terms  of  R  matrices. 


(A- 19) 


(A-20) 


Remember  that  the  C  matrix  transforms  the  vector  expressed  in  column  vector  form  such  as 

r  =  (z|  fr°m  one  frarne  t0  another,  and  the  matrices  representing  subsequent  rotations  are  multiplied 

from  the  left  as  shown  in  Equation  (A-17).  The  R  matrix  transforms  the  vector  form  expressed  in 
row  vector  such  as  rT  (x  y  z)  from  one  frame  to  another,  the  matrices  representing  subsequent  rota¬ 
tions  are  multiplied  from  the  right  as  shown  in  Equation  (A-20). 
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EQUATION  OF  MOTION  FOR  THE  SEMICIRCULAR  CANAL1 


The  semicircular  canal  (SCC)  is  a  rigid  circular  tube  with  a  very  small  uniform  cross-section. 
Each  canal  completes  a  hydrodynamic  circuit  through  the  utricles  (which  it  shares  in  common  with 
the  other  two  canals  on  its  side  of  the  head)  with  negligible  coupling  with  the  other  canals.  The  fluid 
flow  within  the  canal  is  laminar  because  of  the  small,  smooth  bore  of  the  canal.  The  canal  has  high 
viscous  damping  as  indicated  by  the  small  Reynolds  number  of  the  flow.  It  is  common  to  use  lumped 
analysis  and  ignore  the  flow  distribution  within  the  canal.  Thus,  the  fluid  is  considered  to  rotate  as 
a  ring  relative  to  head. 

The  following  analysis  is  based  on  angular  rotation  of  the  canal  in  a  single  plane.  Consider  a 
counterclockwise  rotation  (defined  as  the  positive  direction  of  rotation)  of  the  canal  and  hence  the 
head  as  described  in  Figure  B-l.  The  damping  torque  (positive  when  acting  counterclockwise)  Mu 
on  the  fluid  ring  is  assumed  to  be 

M„  =  -be(F/H)  (B-n 


where  0  (F/H)  is  the  angular  velocity  of  the  fluid  relative  to  the  head  and  b  is  a  positive  proportional¬ 
ity  constant.  Note  that  -b  0  (F/H)  is  clockwise  relative  to  the  head  in  this  case  because  of  inertial 
reaction.  The  elastic  restoring  torque  Me  on  the  fluid  ring  is  assumed  to  be 

Me  =  -k0(F/H)  (B-2) 

where  0  (F/H)  is  the  angular  displacement  of  the  fluid  relative  to  the  head  and  k  is  a  positive  propor¬ 
tionality  constant. 


1.  Reprint  from  Bibliography  17  by  Kee  Soon  Chun. 
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Figure  B-l.  The  Semicircular  Canal  (diagrammatic) 
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By  the  rotational  form  of  Newton’s  second  law, 

J0(F/I)  =  Md  +  Me  =  -b  0  (F/H)  -  k0(F/H)  ^ 

in  which  J  is  the  fluid  ring  moment  of  inertia  and  0(F /I)  is  the  angular  acceleration  of  the  fluid  rela¬ 
tive  to  inertial  space.  By  denoting  the  angular  acceleration  of  the  head  relative  to  space  by  0(H/I), 
it  follows  that: 


0(F/I)  =  0(H/I)  +  0(F/H) 


(B-4) 


Assuming  that  the  cupular  deflection  relative  to  the  head  0(C/H)  is  proportional  to  the  angular 
displacement  of  the  head  relative  to  fluid  0(H/F) , 


0(C/H)  =  a0(H/F)  =  -a0(F/H)  (B-5) 

where  a  is  a  positive  constant,  which  reflects  the  fact  that  the  area  of  the  canal  is  not  equal  to  the  area 
of  the  ampulla.  Substituting  Equation  (B-4)  and  Equation  (B-5)  into  Equation  (B-3) 

J0(c/H)  +  b  e  (C/H)  +  k0(C/H)  =  aJ0(H/l)  (B-6) 


Replacing  cupular  deflection  0(C/H)  by  <|)  and  0  (H/I)  by  H  in  Equation  (B-6), 

J<J)  +  b<{)  +  kcj>  =  aJH  (B-7) 

It  is  the  same  form  as  the  differential  equation  for  a  torsion  pendulum  which  describes  the  tor¬ 
sional  vibrations  of  an  elastic  shaft  with  a  circular  rotor  rigidly  attached  to  it.  The  canal  system  is 
probably  ten  times  more  than  critically  damped. 


Since  the  system  represented  by  Equation  (B-7)  is  very  overdamped  (b2  »  Jk)  and  therefore 
b  J 

t>>i  ,  the  transfer  function  (in  Laplace  transform  notation)  of  the  SCC  may  be  approximated 
by: 

<t>(s)  _  aTlT2 

H(s)  “  (T,s  +  1)  (T2s  +  1)  ^  8) 
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where 


Ti 


T2  “ 


J 

b 


(B-9) 


For  humans,  Tj  is  thought  to  be  about  10  ~  16  sec  and  T2  about  0.003  sec.  For  the  frequency 
range  below  5  Hz  corresponding  to  normal  daily  activities,  where  T2  s  <  1 ,  the  transfer  function  may 
be  further  approximated  as 


4>(s)  _  aT1T2 
H(s)  “  Tjs  +  1 


(B-10) 


Further,  in  the  range  above  about  0.05  Hz,  where  sTi  >  1,  which  includes  most  normal  head 
movements. 


<Ks) 

H(s) 


-  ("2)  1 


which  is  a  pure  integrator  with  gain  (aT2). 


It  follows  that  in  this  frequency  range, 


<t>  (t)  « 


* 

- 


H(t)dt  oc  h  (t) 


(B-ll) 


(B-12) 


Thus,  the  output  of  the  SCC  is  proportional  to  head  velocity  over  the  range  of  natural  head 
movements,  the  role  of  SCC  being  that  of  an  integrating  accelerometer  or  velocity  transducer. 
Indeed,  the  experimentally  measured  firing  rate  modulation  of  primary  vestibular  afferents  is 
proportional  to  head  velocity  over  this  frequency  range. 

The  gain  constants  of  the  internal  neural  signal  processing  elements  are,  for  practical  purposes, 
all  indeterminable  because  the  signals  all  consist  of  the  firing  rates  of  large  populations  of  neurons, 
only  a  few  of  which  can  be  observed  at  any  one  time.  Only  the  final  overall  gain  is  important  and 
interior  gains  may  be  adjusted  arbitrarily  so  long  as  the  total  gain  is  kept  correct.  Thus,  the  gain  of 
the  transfer  function  Equation  (B-10)  of  the  SCC  may  be  arbitrarily  adjusted  so  that  Equation  (B-10) 
will  behave  like  a  pure  integrator  with  a  unity  scale  factor  at  midband  frequencies.  That  is,  Ti  is 
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relabeled  the  cupula  long-time  constant  Tc  and  (aT  1T2)  is  replaced  by  Tc.  It  follows  from  Equa¬ 
tion  (B-10)  that: 


4>(s)  _  Tc 
H(s)  ”  sTc  +  1 


(B-13) 


which  reduces  to  1/s  for  sTc  >  1.  The  break  frequency  for  l/(sTc  +  1)  is  about  0.016  Hz  with  Tc  = 
10  sec  (for  humans).  Thus,  for  all  frequencies  of  normal  head  rotation  above  0.016  Hz,  the  canals 
integrate  head  acceleration  and  produce  a  signal  proportional  to  head  velocity. 
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APPENDIX  C 

EQUATION  OF  MOTION  FOR  THREE-DIMENSIONAL  EYE  ROTATIONS 
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EQUATION  OF  MOTION  FOR  THREE-DIMENSIONAL 
EYE  ROTATIONS 


We  assume  the  Eye  globe  is  a  rigid  sphere  with  a  uniform  density. 

We  use  the  standard  Eye  frame  (Frame  E)  with  x-axis  forward,  y-axis  leftward,  and  z-axis 
upward  with  the  origin  at  the  center  (which  coincides  with  the  center  of  mass)  of  the  sphere. 

We  can  easily  see  that  each  axis  is  the  axis  of  symmetry  in  the  sense  that  the  rotation  of  the 
sphere  about  respective  axis  alter  neither  shape  nor  density  distribution  of  the  sphere. 

This  symmetry  considerably  simplifies  the  derivation  of  the  equation  of  motion.  Although  the 
equation  is  well-established,  it  is  derived  here  from  scratch  to  familiarize  readers  with  vector  opera¬ 
tions  and  with  the  application  of  the  Coriolis  Law. 


By  definition,  the  angular  momentum  (moment  of  momentum)  He  about  C  (the  center  of 
mass)  of  a  point  mass  m^  located  at  point  k  with  distance  Rck  from  point  C  is  (summing  for  all  k’s): 


Hr 


-  ?(Rct 


mkPi 


RCk) 


(C-l) 


d.  j 

where  PjRCk  =  ^RCk  ,  with  -jj  ()  performed  in  inertial  space.  By  definition,  the  torque  Me 
about  the  center  of  mass  C  is 


Mr 


=  ?RCk  x  F 


(C-2) 


where  Fk  is  the  external  force  exerted  on  the  point  k. 

Differentiating  Equation  (C-l)  with  respect  to  time  in  inertial  space  (using  chain  rule): 


>lHc  =  ?  (Pi 


RCk  x  nik  Pj  RCk  +  RCk  x  mK 


Pi  Rck) 


?  (RCk  x  mk  Pj  RCk) 


(C-3) 


because 

PiRCk  x  mi^RQ,  =  mK  (PiRCk  x  PiRck)  = 


(C-4) 
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since  the  mass  mj-  is  a  scalar. 

Now  (by  vector  addition), 

Rfc  =  Ric  +  Rck  (C-5) 


where  I  is  the  origin  of  inertial  frame,  as  shown  in  Figure  C-l  below. 


From  Equation  (C-5): 

Rck  =  Rlk  “  Ric 


(C-6) 
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Substituting  Equation  (C-6)  into  Equation  (C-3): 

PiRc  =  ^  [RCk  x  mkPI  (Rik  _  Ric) 


V"i  V  > 

=  k  (Rck  X  mK  P*  Rlk)  ~  k  (R< 

=  ?(Rck><Fk) 


rr-CkXI»kI>IRIc) 


mk  RCk  I  X  p2  ^CK 


(C-7) 


because  mkPj  RIk  =  Fk  by  Newton’s  second  law. 


Now,  mkRCk  =  0  by  definition  of  the  center  of  mass  relative  to  the  point  C,  which  is  also 

JC 


^  ^  mkRCk  =  0 
the  origin  of  Frame  E. 

It  follows  from  Equation  (C-7)  and  Equation  (C-2): 
PlHC  =  5  (RCk  X  Rk)  =  Mc 


(C-8) 


In  other  words,  the  applied  torque  about  the  center  of  mass  is  equal  to  the  time  rate  of  change 
of  the  angular  momentum  around  the  center  of  mass  with  respect  to  inertial  space. 

To  determine  He  (the  angular  momentum  about  the  center  of  mass)  in  terms  of  the  coordinates 
Xk,  yk  and  Zk  of  Rk,  we  start  with  Equation  (C-l),  which  is  repeated  below: 


He  ~  ^  (mk  Rck  x  Pi  Rck) 


(C-9) 


According  to  the  Coriolis  Law,  treating  the  Head  frame  as  a  fixed  or  reference  frame  and  Eye 
frame  as  a  rotating  frame,  we  have: 


PH  RCk  ~  PE  RCk  +  WHE  X  RCk 


(C-10) 
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where 


ph  RCk  =  ^RCk  ^  Frame  H 


PERCk  =  ^RCk  in  Frame  E 


w 


HE 


angular  velocity  of  Frame  E  relative  to  Frame  H. 


In  Equation  (C-10),  Pe  Rck  =  0  because  Rck  in  eyeball  is  fixed.  It  follows  from  Equa¬ 
tion  (C-10): 


PH  RCk  “  WHE  X  RCk 


Substituting  Equation  (C-ll)  into  Equation  (C-9): 

Hc=?  [mk  RCk  X  (WHE  X  RCk )] 

where  we  replaced  pi  RCk  by  PH  RCk  with  good  approximation. 

From  a  standard  mathematical  table,  we  recall,  for  vectors  Vi,  V2  and  V3: 
Vt  x  (V2  x  V3)  =  (Vj  •  V3)V2  -  (Vj  •  V2)V3 


(C-ll) 


(C-12) 


(C-13) 


Applying  Equation  (C-13)  to  Equation  (C-12),  identifying  Vi  with  m^  Rck,  V2  with  wjjg,  and 
V3  with  Rck^ 


HC  -  (mk  RCk  '  RCk)  WHE  “  5  (mk  RCk  *  wHe)  RCk 


(C-14) 


where 


RCk  “  xk^  +  ykJ  +  zkR 


yk 

ZK 


wHE  =  wxI  +  WyJ  +  WZK  = 


rw? 

w^ 

w- 


dE 

KCk 


=  W 


HE 


(C-15) 


(C-16) 


where  the  superscript  E  implied  that  the  components  are  expressed  in  Frame  E,  or  in  Eye  frame. 
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If  coordinate  axes  are  symmetric  axes  such  as  the  axes  of  Frame  E  fixed  to  the  Eye  globe,  the 
expression  for  jjE  (with  components  expressed  in  Eye  frame)  expressed  in  terms  of  Equation  (C-15) 

and  Equation  (C-16)  becomes  much  simpler. 

Substituting  Equation  (C-15)  and  Equation  (C-16)  in  Equation  (C-14),  we  have,  after  some 
algebra  cancels  out  terms  involving  x^yk,  ykZk,  and  XkZk: 


Hx 

_Xmk(yk +  zk) 

wx 

Hc  = 

Hy 

Hz 

=  Xmk(Xk  +  zk)t 
Xmk(xk +  y?) 

Wy 

wz 

Jx  wx 
Jy  W y 

Jz  wz 


(C-17) 


where  the  meaning  of  Jx,  Jy,  and  Jz  should  be  obvious  from  the  first  equation  of  Equation  (C-17). 
It  follows: 


(C-18) 


or 


He 


r  -i 

E 

Hx 

i 

X 

o 

o 

wx 

Hy 

= 

0  Jy  0 

Wy 

Hz 

0  .0  Jz 

wz 

(C-19) 


where  J  represents  the  coefficient  matrix  on  the  right  side  of  Equation  (C-18).  We  use  J  instead  of 
I  to  avoid  confusion  with  Identity  matrix. 


Equation  (C-17)  to  Equation  (C-19)  show  that  the  direction  of  the  angular  momentum  is  the 
same  as  that  of  the  angular  velocity  for  Eye  frame.  However,  in  general,  the  angular  momentum 
He  and  the  angular  velocity  w  of  a  rigid  body  are  not  in  the  same  direction.  By  definition,  the  three 
axes  of  an  orthogonal  (mutually  perpendicular)  frame  with  this  property  are  called  principal  axes 
of  inertia  or,  briefly,  principal  axes  of  the  body. 
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The  symmetric  axes,  such  as  the  axes  of  Eyeframe,  are  always  principal  axes.  However,  not 
all  principal  axes  are  symmetric  axes.  That  is,  there  are  principal  axes  that  are  not  symmetric. 

Now,  returning  to  Equation  (C-8),  which  is  repeated  below: 

Pi  Hc  =  Mc  (C-20) 

In  our  study  of  eye  rotations,  we  may  regard  the  Head  frame  as  inertial  frame,  with  good 
approximation.  So  we  may  write  Equation  (C-20)  as 

PH  Hc  =  Mc  (C-21) 

Applying  the  theorem  of  Coriolis  between  Head  frame  and  Eye  frame  for  PhHc  in  Equa¬ 
tion  (C-21): 


Ph^c  —  Pe^c  +  whe  x  He  —  Mc 


(C-22) 


In  other  words,  the  applied  torque  about  the  center  of  mass  of  eye  produces  and  is  equal  to  the 
rate  of  change  of  angular  momentum  with  respect  to  the  Eye  frame  plus  the  cross  product  of  the  an¬ 
gular  velocity  of  the  eye  relative  to  the  head  and  the  angular  momentum  of  the  head  about  the  center 
of  mass. 


Using  notation  given  in  Equation  (C-16)  and  Equation  (C-17): 


WHE  x  Hc  - 


I  J  K 

W x  W  y 

hxhyhz 


I 

J 

K 

= 

wx 

wY 

wz 

Jxwx 

JyWy 

Jzwz 

fjzwYwZ 

J  yW  YWZ 

= 

1 

N 

>—> 

Jzwxwz 

1  JyW^W  y  — 

JxWyWy 

(C-23) 
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and 


PeHc  = 


PE 


Hy 

Hz 


Jxwx 

JYWY 

Jzwz 

Jxex 

Jx  ®x 

JyQy 

= 

Jy  0y 

Jz6z 

Jz  ®z 

(C- 24) 


(C-25) 


Substituting  Equation  (C-23)  and  Equation  (C-24)  into  Equation  (C-22)  and  using 
Mc  =  Mq  =  [Mx  My  M2]  : 

Jx  ^  wx  +  (Jz  ~  Jy)  wywz  =  Mx 


h  (jt  wy  +  (h  h)  wxwz  ~  My 


(C-26) 


Jz  ^  wz  +  (Jy  ”  *x)  wxwy  - 


Equation  (C-22)  and  Equation  (C-26)  apply  to  the  body-fixed  (eye-fixed  in  our  case),  princi¬ 
pal  axes  of  a  rigid  body  (Eye  frame  axes  that  are  symmetric  in  our  case)  with  origin  at  the  center 
of  mass.  In  general,  these  equations  represent  highly  nonlinear  differential  equations  that  are  diffi¬ 
cult  to  solve.  The  solution,  when  obtained,  is  the  angular  velocity  of  Eye  globe  with  respect  to  Head 
expressed  in  Eye  frame,  and  does  not  give  directly  the  motion  relative  to  Head  axes. 
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For  a  special  case  of  only  an  axis  rotation,  such  as  x-axis  rotation,  wx  =  wY 
from  the  first  equation  of  (C-26), 

j 

JX  Jt  wx  =  Mx 
d0y 

Denoting,  wx  =  dy  -p-  =0  ,  we  get  from  Equation  (C-27)  that 

dt  x 

Jx  0X  =  Mx 

which  is  commonly  used  in  the  one-dimensional  analysis. 


0.  It  follows 

(C-27) 

(C-28) 
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